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INTRODUCTION 


Fracture Problems and Fracture Criteria; A Review 
Character lzln(i the dleplacement, atrees and strain fields associated 
with stationary and propagating cracks in solids characterized by various 
idealized constitutive relations is one of the most Important areas of study 
in fracture mechanics. The Importance of these studies Is not only that we can 
predict the stress or displacement fields In a cracked body, but also that 
knowing the nature of such fields we can possibly correlate observed fracture 
behavior with some aspect of these fields and thus arrive at valid fracture 
criteria. 

Criteria for Crack Growth Initiation 

The two macroscopic aspects of fracture for which correlations are commonly 
sought are the initiation of crack growth and the rate of crack growth. The 
most notable initiation correlations are with the elastic stress Intensity 
factor, K^, for the elastic (and/or small scale yielding) case [1,2] and with 
the Jj^-integral for cases In which plasticity may not be limited to the ctack- 
tlp [3,4]. The conditions under which these correlations are independent of* 
geometry are discussed In the cited references. The critical values of and 
for a given material are denoted and respectively. It is implied 

by the use of the subscripts "I" and "1" that these criteria are for the crack 



opening mode (l.c.» mode 1). While similar criteria might be expected for the 
sliding und tearing modes (l.e.» modes II end III), the exper^aental data for 
such studies is lacking. 

Both and have been shown to be crack-tip field parameters and both 
can bo evaluated experimentally through energy considerations. The linear elas- 
tic, mode I civ'.ck-tip field determined in [5] shows that the asymptotic crack- 
tip fields arc proportional to K^. Similarly, the asymptotic, mode I, crack- 
tip fields for power-law deformation theory plasticity have been shown [C. 7] to 
depend upon the single parameters In the case of small scale yielding, 

Is easily related to the energy release rate, [8], which la a measure of the 
potential energy decrease due to an increase in crack length. The quantity 
has a similar potential energy interpretation in the case of deformation theory 
plasticity, and becomes identical to for small scale yielding of a stationary 
crack. 

In the foregoing discussion, the time dependence of the material's response 
and of the applied loading is assumed to be negligible. For creep crack growth 
those assumptions arc no longer valid. We now consider crack growth initiation 
in materials which exhibit creep behavior. While a significant number of creep 
fracture experiments have been reported in the literature, it appears that the 
primary interest has been to find a creep crack growth rate criterion as opposed 
to an initiation criterion. As a result of this emphasis, many investigators 
use notched specimens rather than precracked specimens and many do not report 
data wiiich could be useful in addressing the question of init iation. At 
present there seems to be some indication [9,10] that when precracked specimens 
are used, the time required for creep crack growtli initiation is negligible when 

^ The deformation theory of plasticity precludes elastic unloading from an elastic- 
pla Stic state and thus Is mathematically equivalent to nonlinear elasticity. The 
crack-tip fields associated with power-law deformation plasticity are commonly 
referred to as HRR f'elds after the authors *f references [6] and [7). 
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compared to the life of tho spucimen. It should be understood that this is not 
a universally acknowledged concluoion [11] and that further atudy la Indicated. 

As noted previously, the second macroscopic aspect of fracture for which 
correlations with crack-tip fields are coimnonly sought Is crack growth or propa- 
gation. The following summarizes the development of criteria relating to this 
aspect . 

Criteria f:u‘ Crack Propaeatlon 

Slow crack growth occurring under constant load Implies that the material 
response Is time dependent and Is generally classified as creep crack growth. 

If the material's time dependent nature Is negligible under the subject con- 
ditions, then It is assumed that cra';l'. growth requires an Increase in applied 
load. This latter case is typical of situations in which small scale yield 
conditions are not met and for which has been found to correlate with crack 
growth initiation. The primary Interest In this quasl-statlc mode of crack growth 
is that for some materials and geometries, the increase in load carrying capacity 
of the structure during quasl-statlc crack growth Is significant. This implies 
that design procedures can be developed to take advantage of this added margin 
of safety. To justify such a procedure, however, there must be some dependable 
means of predicting the crack growth versus load behavior as well as predicting 
at what load the crack becomes unstable (l.e., is no longer quasl-statlc). As 
noted, is generally accepted as a valid initiation criterion for this problem. 
For prediction of the subsequent growth, however, there are at least two proposed 
criteria which appear to provide reasonable correlations with experimental data. 
The first growth criterion can be stated as , where 

JirCas) is assumed to be a material property which depends on the amount of 
crack growth, a [12,13]. The subscript "R" denotes that this quantity 
characterizes the material's resistance to cracking. While strong theoretical 
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argumenta can bo glvan aa to why this criterion should not be valid (except 
possibly for very United amounts of crack growth (14]) » it has been denon- 
atrated that reasonable predictions can result from the use of this criterion 
for at least some classes of problems [12,15]. 

Based on the theoretical objections to the use of except for United 
crack growth, a second criterion which is based on the crack-tip deformation 
hjs been proposed [16,17]. This criterion results from finite element simu- 
lations of quasi-static crack growth experiments which indicate that the crack- 
tip opening angle, CTOA (defined by the first finite element behind the crack- 
tip) becomes constant during crack growth. Whereas the CTOA, so defined, is 
clearly a mesh dependent quantity, the concept of crack-tip deformations becoming 
constant with crack growth is physically meaningful. The procedure for applying 
this criterion in finite element based predictions of crack growth behavior is 
to use for Initiation and for crack growth prediction until the computed 
CTOA has become constant with crack growth. Continued growth is then governed 
by this constant value of CTOA. Alternatively, a predetermined CTOA resistance 
curve can be used throughout growth. Crack growth instability is assumed to oc- 
cur (for either CTOA or J, as the criterion) when further increase in crack 
length results in the criterion for growth being exceeded without further increase 
in applied loading. The J._ and CTOA criteria appear to provide reasonable 
correlation of ductile slow crack growth behavior for a variety of materials, 
geometries and load conditions [15,18,19]. 

Creep crack growth generally becomes a concern when components are operated 
at elevated temperatures. Whereas quasi-static crack growth can be on the order 
of nun/scc, typical creep crack growth rates arc on the order of pm/sec. Compared 
to clastic-plastic quasi-static crack growth, the problem of creep crack initiation 
and growth is a relatively new area of study. 
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Numerous expsrlmsntsl studies have been undertskcn with the purpose of 
finding e perenojter which correlates with creep crack propagation rate. (Sea, 
for example, the review article [20] and [21-24].) Most of these investigations 
consider as candidate parameters, K^, some form of net section (or reference) 
stress, and in more recent atudlos C*. The C* parameter is the steady-atate 
creep analogue of (in the sense of [25]) in that the definition of C* is the 

same as that for except displacements and strains are replaced by their res- 
pective rates [26]. 

It is illustrated in Fig. l.I that the above three parameters can be expected 
to correlate three distinctly different creep crack growth situations. In Fig. 
l.Ia, a crack and its associated ligament are shown for a material and geometry 
which results in negligible creep strains everywhere except in the vicinity of 
the crack-tip. This condition is analogous to that of small acale yielding in 
elastic-plastic fracture. Fig. 1.1b represents a situation in which C* might 
be considered an appropriate parameter. This situation is characterized (1) by 
the body being essentially at steady-state creep conditions (which implies very 
slow crack propagation) and (11) by the creep-damage process-zone being local 
to, and therefore controlled by, the crack-tip field. Fig. 1.1c Illustrates 
the type of situation for which net section stress might be expected to control 
crack growth. In this case, the main feature is the widespread creep damage 
zone. 

It is seen from Fig. 1.1 that intermediate situations can occur. For 
example, suppose a particular material and geometry results in a crack propagation 
rate such that elastic strain rates are not negligible compared to creep strain 
rates (l.e., nonsteady creep) and at the same time, creep strains are no longer 
localized to the crack-tip region. VHiile neither or C* could be valid 
parameters for this case, it apoears reasonable to expect that crack growth 
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Fig. 1.1 Conditions for which creep crack growth 
parameters are expected to be valid 
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rat« it ■till d«t«raln«d by tl)« local ccack-tlp field alnco the craap daaaga 
procaas sona la atlll aaauaMd to ba local to the crack-tip. 

A paranatar which la apparmtly capable of apannlng tha t«P batwaan 

controlled growth and C* controlled growth haa baan Introducad [27]. Thla 

paramatar la rafarrad to aa (AT)^ and la defined by a path* Independent vector 

integral. A detailed dlacuaalon of a ganarallaad C* (l.e.. C*) and 

la given in Section 11 of thla work with a principal reault being that the 

energy relatlonahlp coaanonly uaad for axparlaental meaaurement of C* doaa 

not apply to C* but rather appllaa to tha (A^) parameter. Thla maana chat 

the experimental reaulta are actually ahowlng a correlation with (AT) rather 

c 

than with C*. baaed on the theoretical validity of (AT)^ aa a crack-tip field 
parameter for nonateady aa well aa ateady-atate creep and baaed on the mounting 
.xp.rla.nt.l «W.n« th.t cr.ck prop.g.tlon r.c. crr.Lt.. «.U with (4T)^. 

It aeema the creep crack growth rate problem la cloae to having a aolutlon. 

Motivation for tha Preaent Work 

In the foliowlngp we review previoua atudlea to the extent required to 
place the preaent atudy in perapectlve and briefly Introduce the preaant work. 

The nonlinear nature of creep conatltutive relations precludea analytical aolu- 
tlona for cither stationary or propagating cracke In a creeping material. For 
stationary cracks in a power*^)aw creep material, however, it la known that 
the HRR fielde are preaant In the vicinity of the crack-tip 126]. (Since the 
singularity In creep strain rates is greater than tluit In the atreaaea, and thus 
elastic strains. It follows tliat the HRR field exists at the crack-tip during 
nonstcady as well as sti^ady-state creep.) For propagating cracks , It appears that 
the HRR fields no longer exist at the crack-tip, but that analytical tools exist to determine 
the fields which do exist [28, ]9]. Wliile knowledge of the crack-tip field is 
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valuftblc, th« lolutlon of boundary value problMu* auat dopond on nuaarlf.al 
mathods. Tho finite aloMnt aathod, In particular, ahowc pronioa for aolving 
creep crack growth probleaa. 

Only a few atudlea on finite element modeling of creep crack growth have 
been reported. The earlieat ie apparently that of Ohtani and Nakamura (30). 
Thie Btudy eimulated crack growth with .« noda-releaae technique and aeeumed a 
critical crack-tip plaetic etrain criterion for creep crack growth. The rate 
conetitutive law contained an elaetic term and a creep term baeed on the gen- 
eralization of the uniaxial Norton power law. 

Hinnercha [31] uaea the Bodner-Partom conetitutive law [32] and a noda- 
releaee technique for modeling crack growth. In thie work, eeveral candidate 
criteria are examined by aimulating crack growth experimente. Due to the ap- 
parently limited crack growth (<0.5 mm), the abort teat duratlone (one hour) 
and the lack of crack growth meaeurement data (which requirea the developaent 
of a ao-called hybrid experimental-numerical procedure to eetlMte the crack 
growth hiatory), it aeema the general applicability of the concluaiona from 
thia otudy are queatlonable. It alao aeema likely that the methoda for evalu- 
ating In thia atudy are incorrect^ and thua the concluaiona concerning C* 
ahould be regarded accordingly. 

Ehlera and Riedel [33] have conducted a finite element analyaia of a 
atatlonary crack in a compact apeclmen. The primary emphaaia in thia atudy 
la on the nature of the crack-tip field during the tranaitlon from the initial 
elaetic field to the ateady-atate creep field. 


While the detaila of the numerical procedurea for evaluating C* are not given 
in [31], It appeara that the W* term of C* (aee Section II) ia incorrectly in- 
terpreted aa a hiatory dependent quantity aa oppoaed to a quantity dependent 
aolely on the atcady-atate atreaa and atrain rate. 
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OF POOR QUALITY 

Tlu^ mri'iiBih of tlu^ lIRR field <lurLni5 Uh' tr<mBltlon period la determined 
tltrou;^!) fitting the neat’ tip equivalent Btress fluid. The calculatlonis u»o 
elliht-nuded laoparmnetrlc elumenlB with quarter-point clemcntB belnft used at the 
crack-tip so as to have an r “ “ strain slnRtilnrlti . Creep crack growth and 
creep crack growth criteria are not considered In this study. 

The finite clement equations for the creep crack growth model being used in 
the present study are derived from the principle of virtual work in Section III. 
Section V presents the rcBults of several analyses involving both stationary cracks 
and propagating cracks. The creep crack growth simulation Is via a mesh shifting/ 
remeshing procedure. Calculations are made using ♦'he quarter-point element tech- 
nique as well as with a specially developed (Section tV) compatible element which 
Incorporates the HRR, ^ strain singularity. 

Ail important aspect of the current work is the study of the (AT) parameter. 

In particular, the meaning of its relationship to C*, and its calcula- 

tion within the context of finite element analysts arc explored in depth. 

A scrlcc of crack propagation calculations are combined with analytical and 
experimental renults in Section V to show that creep crack growth in 304 stain- 
less sLi'cl at 6!)0^C occurs under essentially steady-state creep conditions. This 
implies that the crack grov.i.h rate for a given crack length and load can be 
detcnnlned from a steady-state t. oep solution which does not depend on the 
previous load and crack growth hl8tor,.eH. This observation implies that 
simple crack growth prediction methodologies irniy be developed. 
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SECTION II 


DERIVATION OF THE (AT>^, J AMD C* INTEGRALS 


Pr«llttlnarl»i 

Wt ahall consider problems which exhibit the following constitutive be- 
haviors 


IJ " *1J ffjj 


( 2 . 1 ) 


We denote the cartealan coordinates of the undeformed body as x^. Defining 
as ths rate of displacement (or velocity) of a isaterlal particle from the cur- 


rent configurer Ion » then c.. is the symmetric part of the rate of displacement 
gradient “ij * gradient operator V^ Is with res- 


pact to the current coordinates y^ where it Is understood that ^ 

^IJkt tenaor of Instantaneous slastlc moduli. We let Oj^^^ denote the 


Ij 


corotatlonal rate (or "Zaremba-Jaumann rate") of the Klrchhoff stress o 

3y_ 

where Is related to the Cauchy stress by - JT^j(J«det[-5jp]) . 

equivalent Klrchhoff atress o Is related to the devlatorlc Klrchhoff stress 
’ eq 


(- 0 ,, - 1/3 0 .. 6 ,,) by o_ - (3/2)(o' o' ) 


1/2 


ij ' •'ij " '^kk''lj' ' eq 

those of the familiar Norton's law 


U lj‘ 


The parameters Y and n are 


€ - y(o )” 

eq eq 


were 




1/2 
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W« will us« tht notationt (.) danotas a aacond ordar tanaor; (_) inpllai 
a vaetor; l • » . £ inpliaa a^ - S l»plU§ 

A:B • Alao nota that £^BtC lapllaa and . | inpliaa 



A Conaarvation Law for Finita Ilaatic and 
Wonataadv Craan Hatarial Bahavior 

rha diacovary of conaarvation lawa and tha poaaibiilty of daclving path- 
indapandant intagrala from thaaa lawa ara not particularly racant occurrancaa aa 
diacuaaad in [34]. Howavar, tha litaratura in thia araa haa baan rathar piaca- 
aaal and tharafora difficult to aaoinilata. Tha racant work of Atluri [27] haa 
dona much to unify and ganaralina thia aubjact and ia tha baaia for tha following 
praaantation. 

Wa will conaidar a vary ganaral conaarvation law which haa baan givan by 
Atluri, but will limit our diacuaaion of thia law to aatarlala chaxactarisad 
by (2.1). Ua will uaa cartaalan coordlnataa axclualvaly. Nota that by 
apacial aalactlon of matarial conatanta (i.a. , Y”0), (2.1) can ba apeclaliaad 
to alaaticity. Altamativaly, by aaauming that tha atraaaaa ara invariant with 
time, (2.1) can ba apaciallzad to ataady-atata craap bahavior. 

In the following praaantation, the currant configuration (i.a., tha config- 
uration at time t) ia tha rafaranca configuration. There may ba initial atraaaaa 
axlatlng for thia rafaranca configuration. If atraaaaa do cxiat, than they ara 
aaaumad to aatlafy tha linear and angular momentum balance condition (I.a., 
equilibrium) 


7^ . T + p^(f^ - o^) - 0; I ■ 

where p^, ^ and a^ are the current maaa dei\^lty, body force vector and accel- 
eration vector. 
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A conf«rv«tion inttgral r«l«tlon glvtn by Atlurl [27] for • elotod volumo 
(ot tho currant tlaa, t), which la fraa froa alngularltlaa and any othar 
dafacts (which would pracluda tha application of tha dlvarganca thaoram), la: 


0 - 



(7 dW - (7 t):A| - 7 


Ul+At) . Aa] 


( 2 . 2 ) 


- P.(f-a) . Aa)dV + /* liv . (t+At) - tl . AadS 

t - - Js, ~ 

0 


+ /* IL. • (^■♦'At) . (Aa-Aa)dS 

a 

In (2.2), At la tha Incraaantal f Irat-Plola-Klrchhoff (nonaynmatrlc) atraaa 
(At ■ [Ag - Aa . g]/J) whara Ag la tha matarlal Incramant of Klrchhoff atraaa. 
Tbo current maaa denalty la denoted and f^ and a are the body force and 
acceleration vectora at time t+At, reapectlvely. S and S are the port Iona of 
the boundary of upon which preacribed tractions, £, are acting and at which 

prescribed dlaplacement gradients, Ae, exist, respectively. The current outward 
normal to or is i^. The quantity AW, discussed in detail in [27], is the 
incremental stress-working density in time At, and is given by: 

AW ■ x:Ae + ~At^:Ae = x:Ae + AU (2.3) 

where 

AU - ^t^iAe (2.4) 

The validity of (2.2) is readily verified through the two identities [27]: 

7 AW - ^(t:Ac) + ^AU - ^x:e (2.5) 

— t — t - — t -H' - ‘ 

T 

+ 7 Ae: X + ^Ae:At 
— t ' ' “t ' 

and 

^ . [(x+At) . Ae] - [7 . (x+At)] . Ae + 7 Ae:(x+At)^ 


( 2 . 6 ) 


ORIGINAL PAGE IS 
OF POOR 'QUALITY 


the ■•titfftcelon of llnoar nomontus bolaneo In V^t 

^ . IT+Atl + - 0 (2.7) 

and tha aatlafaetlon of tha boundary eondltlona:^ 

n^ . [t+dtl • E on (2.8) 

da « da on S (2.9) 

Nota that Idantlty (2.5) aaaumaa that t ( tha initial atraaa for tha incra- 
mant) la an axpliclt function of Ita pooltlon In V^. Tha axiatanca of dU 
la ahown and dlacuaaed In tha work of Atlurl [..'5]. 

Having tha ralatlon (2.2) It la now poaalbla to apaciallza thla ralatlon 
to finite alaatlc behavior or to ateady>atata creep behavior. However, alnce 
we are primarily Interested In the path- Independent Integrals which can be ob- 
tained from (2.2) we will postpone the specialization till after we have de- 
rived the general path-independent Integral (dT^) . 

Path-Independent Integrals for Fracture Analysis 
The conservation Integral (2.2) Is used [27] to obtain a path-independent 
integral which Is applicable to the analysis of cracks by considering a volume 
such as Illustrated In Fig. 2.1. (Note that a two-dimensional case 
Is Illustrated for simplicity) . The use of the divergence theorem for the 
region depicted In Fig. 2.1 results in (2.2) being rewritten 


+ 



[n dW - n . (t+dt) . de]dS 
t — t 


[(-7 T);de -p (f-a) . de]dV 
”t* ' t 


( 2 . 10 ) 


The validity of (2.2) does not require S^'fSg^dV where dV denotes the surface 
bounding V^. Therefore, dV need not coincide with the boundary of the body 
under consideration. 
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+ f n^AWdS + J n^AWdS - J t . "*dS 

’ J SL^ • (T+At) . Aids 
“ J IHj^W - . (t+At) . AeJdS = (AT)* 


( 2 . 10 ) 


In writing (2.10) it has baan axaumad that ^12^^45* inpliea that 

T234 coincida with any axtarior boundariea. Thia haa baan aaaumad 

purely for conveniance of notation. We have also used the notation F •• F,,.. 

e 165 

Noting that (2.10) contains two equalities, it c.«n bn verified by inspection 
that (AT)* depends on e (or more generally F ) but that it does not depend 
on the selection of ^' 234 * this sense (AT)* is path-independent (i.e., 
independent of the selected fsr-field path). Following the reasoning of Atlurl 
[27], we define (AT) as the limit of (AT)* as € goes to zero.^ 


(4X) - 

^ -''c €-*0 


/ [n AW - n . (l+At) . AeJdS 
-t -t - 


( 2 . 11 ) 


' -^234 


[n AW - n , (t+At) . Ae]dS 


lio) f 


+ f nj,AWdS + r n^AWdS - T t . AedS 
*'^2 *^^5 ® 

- /* n . (r+At) . Aeds! 

Js ( 


The existence of the limit is shown in Appendix A. 
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Pig. 2.1 Contours for applying the conservation lav/ to 
a two-dimensional, cracked body 
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By defining (AT)^ ae the limit of (Al)* ■■ « io«* to tero» It le seen that 

(&T)^ la entirely determined by the aeymptotlc near-tip flalda. It will be 

shown later that the converse to this atatement la also valid when the near- 

tlp fields are the HRR flelde. That la, It will be shown that (dT) entirely 

c 

determines the asymptotic near-tlp fields. 

Often It happens that only the flrat component of the vector quantity (^T)^ 
la of Interest. We will write the first component of (dT) as (dT. ) . dlao, 

C X c 

the quantities (T^)^ and will often bo used in place of (dT)^ and 

These quantities are related by 


c 


U 

df*-0 dt 


( 2 . 12 ) 


However, in the presentation of numerical evaluations of (dT.) we use (t. ) 

X C X c 

as a convenient, approximate notation for (AT.) /dt. 

i c 


We now consider two special cases of (2.11). For symmetrical deformation 
about the x^ axis and cracks oriented along the x^ axis with traction free 
crock surfaces, no body forces and negligible inertial effects, the first 
component of (dT)^ is 


(AT^) 


c 


Lt 

£-►0 



[njdW - nj(Tj^+Atj j^)Ae^j^]dS 


(2.13) 


- I In^AW - 

^234 

- f dV 

J 3yt 

Note that the limit of the volume Integral has been written In Its explicit* 
form as a result of the arguments for the existence of this limit, given In 
Appendix A. If, in addition to the above conditions, the strains are in- 
finitesimal and the deformations small, then there is no need to differentiate 
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between and l>*comee Identical to *nd we have (recalling 

Ac^j r. DAu^/3yj): 

Lt r 

- Jo J, IV" - «•“> 

« ^ 

3Au 

|nj4W - -j^ldS 

234 ^ 

- 4 

t 

The replacement of Ae^ by Ae^ in the volume Integral is made posaible by 
the symmetry of and thus 3x^j/dXj^. 

Physical Interpretation of (AjT)^ 

It has been shown by Atluri [27], that the vector (A^ has the following 

c 

physical meaning. Let two bodies with non-propagating cracks be identical 
except for the second body having an additional, arbitrarily directed, in- 
finitesimal Increment in crack length characterized by the vector It is 

assumed that both bodies experience Identical load histories. Define total 

2 

potential energy increments corresponding to the time increment At as 

AE^ - Atj>^ + Afl^ + AK^ (2.15a) 

AE^ - AiJ;^ + A«2 + AK2 (2.15b) 

for the first and second bodies, respectively. In (2.15), -Ai]> is the incre- 
mental work of externa] forces, LU is the incremental stress-work and AK is 
the increment in the kinetic energy. (It should be noted that AS] includes the 
Inelastically dissipated energy.) Then 

^Atluri [27] has shown that the 1/r singularity in kinetic energy, which is as- 
sociated with dynamically propagating cracks, changes the interpretation. 

2 

Note that sign convention for AE^ and AC is opposite to [27] so as to conform 
to conventional usage. 
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(AT^)^bdc^ - -(AE^ - AE^ 


(2.16) 


whore b iu the length of the crack front. 

If one is only interested in self-similar crack extension in the -direct ion, 


then dc^ • dc^ 


0 and 




^^2 - 
bdc. 


(2.17) 


Therefore, (AT)^ is related to the Incremental potential energy difference 
between two bodies which are identical except for an incremental crack length 
difference 

Finite Elasticity and £ 

As noted previously, the constitutive law (2.1) can be specialized to 
elastic behavior by choosing y to be zero. Therefore, (AT^)^ as defined by 
(2.11) is a valid crack-tip characterizing parameter for general nonlinear 
elasticity with finite strains, large deformation, body forces and inertial 
effects. Howerver, the basic premise of elastic behavior is that the con- 
stitutive relations arc independent of the histories of deformation and stress. 
This means that the constitutive relations can be derived from a potential. 

For instance, a potential, U, exists for the first Plola-Klrchhof f stress, 
such that 


T 3U 
^ " 3e 


(2.18) 


In the following, we consider the reference configuration to be the stress-free, 
undeformed configuration at t"0, and therefore drop the subscript t for con- 
venience. As a result of the existence of the relation (2.18), it is possible to 
state two identities which are analogous to those of (2.5) and (2.6) for the case 
T=0. 


3U mn .. . T 

VU » -r - Ve;t 

— 3e 3x, — >• ' 

mn 1 


(2.19) 
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1 • It • ?1 ■ Z • t*f + Zf‘S (2.20) 

Similarly, the linear momentum balance (i.e., equilibrium) condition le now 

1 • 5 P<l. - •) • 0 (2.21) 

end the boundary conditions 

n . E • 1 on (2.22) 


s • e on S 


(2.23) 


Noting the similarity of equations (2.16) through (2.23) with (2.4) through 
(2.9) it Is easy to arrive at the following conservation law 


(2.24) 


0 - / {vu - V . (t . el - p(f - ^ . eldV 

Jyj 

+ f In • t “ JLl • §<^5 + f n. • t . (e - e)dS 

Js^ “ Js 

t e 

Following the procedure used in deriving (4T)^ from (2.2) we apply the diver- 
gence theorem to (2.24) for the volume V - and take the limit as e goes 
to zero. The path-independent vector quantity resulting from this procedure will 
be called J. 


. <•' f 

c-O J, 


[nU - n . t . e]dS 


(2.25) 


r ( ’• 

■ j [nil - n . t . e]dS t J p(a - O . edV 

^234 

+ / nUdS + f nUdS - f F . edS - f n . t . edS | 

" •'s ” ” A ” ' “ 

1? 45 t e ' . 

In writing (2.22) through (2.25) It la understood that o Is the mass density 
In the reference configuration, £ Is the unit normal In the reference configu- 
ration, £ are arbitrary body forces per unit mass, a Is the absolute material 
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original paqi w 


accuieraclon, and jt arc prascrlbad tractions par unit undafomad area S^. 

Wa now conaldor aaveral apaclal caaaa of (2. 23). If tha problaa balng 
conaldorad involvaa a crack orlantad along tha coordlnata dlractlon and la 
loaded ao that only mode I crack-tip behavior occurs, then la of primary 
Inturcat and we have: 

- X ‘“i" ■ "I'ij’ji’" * L V 

*234 ' ^ € 

- [ '1*11'*® - X 

t e 1 

If In addition the problem Involves infinitesimal strains, small displacements 
and traction-free crack surfaces, we have 


J 


1 



■ ‘'iHj 



(2.27) 


where use has been made of the existence arguments of Appenlx A in taking the 
limit of the volume Integral. 

For elastic behavior and non-propagating cracks, Atlurl [27] shows that 

J has the meaning of energy release rate to a process zone V in the sense that 

c 



where b is the length of the crack front, 

DE Dii) on DK^ 

“Dt" “ “dT "dT 


and 

dt “ E (t + dt) - E (t) 

e c 

For an clasto-dynamically propagat 'ng crack (i.e., singular kinetic energy) 
Atluri [27] concludes that 






DL 

t 

Dt 
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kjkj, IV# 


dt 


(2.29) 


whora L !■ th« Liigrangian L, i X, - auch that 

c c c c c 


DL. 


Dt|« 


DK^ 

“dT 3T 


Dt Dt 

Therefore, haa the meaning of "rate of change of Lagranglan per unit crack 
growth". 

Wo now conaider the apodal caae of ataady-atate creep behavior. 

Stoady-Statc Creep and C* 

It haa been ahown that (AT^) characterizea the crack-tip field for mat- 

c 

erials which exhibit creep behavior auch aa in (2.1). It ia known that under 
certain conditlona of applied loading, the conatitutive relation (2.1) can (after 
long t lilies) result in a steady-atate. Thia ateady-atate ia primarily charcter- 
ized by the time independence of the atreaaes (i.e., AU ■ Specializing 

(2.13) to steady-state conditions, we define the steady-atate value of (AT^j^)^: 

e 

■ f. - X 

‘23A '^t ^ 

Because (2.1) results in a power-law relation at steady-state, which is 
analogous to the power law deformation-theory plasticity (or essentially non- 
linear elasticity), Goldman and Hutchinson [26] have suggested a path- independent 


C* integral, 


c* - j" I 

•'r L 


"i”* " IT, 


dS 


(2.31) 


whore 


W* 


io' 


T ,d^ 

ij ij 


(2.32) 
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The quaatlon of how C* and (<^T ) aro related , ie a natural one. Before 

^ X V V 0 

obtaining an equation rolatlng C| to however, the coneervetion Integral 

(2.2) will bo used to derive a generalUed vector Integral C*. 

In epeciallzing (2.2) to eteady-etate we note that now eireee la a 
function of the strain rate, and that etreei Incremente arc zero. Thus, 

AW'T;Ae. Also wo may write: 



(7 AW - (7 T):AeldV 

“T t - ■» 



T :7 AedV 
- ~t * 


Thus, at stoudy-stato, we may write (2.2) as: 


0 “ f - 2 * • ^T.Ae) - p (f - a) . Ag]dV 


+ y* (n^ . T - jt ] . A£dS r £Lj. • I • - Ag)dS 

t e 

or equivalently, in rate form, 

0 - f (t:(V e) - V . (T.e) - P (f - a) . eldV 

-t t - - 

I llL. • T - t] . edS + f n . T . (e - e)dS 
t c 

Using the symmetry of t we note that: 

t:Vc - T:vf^ (e + e^)l " t:V c 
- -t ' '•~t*-2 ' ~ J ' — t - 

As a result of the incompressibility condition ■ 0) we !iave 

t:V^* - 


and W* ■ 




which leads to 


(2.33) 


(2.34) 


(2.35) 

(2.36) 
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Combining the rosults of thr abov« manlpulmtions w« luivmt 


7 W* 


aw* 

hi 


aw* ^^mn 
ay 

mn *^1 


I'.V J . ;,» * 


Using (2.38) and tha dlvarganca thacram whlla applying (2.33) to 
daflna th# vactor quantity (C*)*: 


(2.37) 


(2.38) 
V^, wa 


f (!LW* - t - •) • * 

Jv ~t - > Jm .V ^ 

‘234 


LdV 


(2.39) 


v^-v 

t € 


+ J n^W*<*S + J n^W*dS • f t . ^dS - J n^ . T . |dS 

•'^^2 •'*'45 •'St •'S, 

- f In W* - n . T . e]dS s (C*)^ 

Jr ~t —t - - 


c 

If we define the limit of (C*)^ aa c-K) to be C*, we have a quantity which 
characterizes the crack-tip field and Is independent of the selection of ^234* 
Restricting our attention to problems involving symmetric deformations about 
the axis and cracks oriented along the axis, with traction-free crack 
faces, no body forces and negligible inertia effects, we find that 


J r* 


(2.40) 


/.. '"i”* - "j 

* A s 


234 

In computing W* it ia convenient to invert (2.1), substitute the re- 
sult into (2.36) and use the following identity to complete the Integration: 


1-fn 

d(c ) n 
-eil_ 


dl 


IJ 


2 ,l+n. f. . 

3 n eq' ij 


(2.41) 
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Ttu r««ult of this Mnlpulatlon in 
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I'fn 

W* . (Y)"(^ ) n 


n y. vH, 

' «q 


or 


W* ■ TIT V(o 
l4f» «q 


(2.42) 

(2.43) 


Rolationohip of for St«ody~Stat« Croop 

Now w« will rolato C| ot (2.a0) tr th« •toady-scato valua of 


Flrat wa rawrita (2.30) in rata form aa: 


\ _ Lt 

cat c-^O 


/p ^"l^ij*ij " "j^ji^il^*^^ 


(2.44) 


234 


Uaing the notation have 




(T.) 

I CIS 


X ‘"i“ - 


il 


)dS 


(2.45) 


234 
3 t 


_ f UL A 

Jy 3yi 


dV 


Noting that: 


“ ■ • 'ij 2 * 'jl> • ^j'tj 

It la seen that W la the rate of atreaa^workina denalty . while W* la luat a 
mathematical potential for t |y Aa a result of Incomprcasibility we can write^ 


(2.46) 


1 


W - 0 • - y(‘J - (“)" ) 

eq eq eq Y eq 


1-Hi 

n 


(2.47) 

t 

as contrasted to W* of (2.42) and (2.43). Comparing the left equalities of 


This result is only valid for steady-state creep and is obtained through the 
substitution of the steady-state specialization of (2.1) into (2.46). 
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(2.40) and (2.44), it U aoan that •"«* CJ ara ralatad by* 






(2.48a) 


" * n+1 €-K) 




P * aq 

€ 


(2.48b) 


Appendix B gives several numerical examples of relation (2.48) for two rather 
extreme values of n. 

It is now clear that £<*> and (4p are not equivalent quantities under 
any condition despite their being derivable from the same conservation law . The 
quantity (4T)^ follows more directly from the ''.onservat ion law and is the more gen- 
eral quantity not only in that it is applicable to nonsteady as well as steady- 
state creep but also in that it is applicable to constitutive laws which are more gen- 
eral than (2.1). The quantity C* relies on the special property of (2.1) which allowa 
the existence of a potential W* for the stresses (t'). Furthermore, since W* does 
not have any physical meaning, whereas W has the meaning of stress-working density, 

it is understandable that (AT) has an energy interpretation whereas C* does not. 

*" c 

It is for this .''eascn that it seems more appropriate to refer to experimental measure- 

dp 

nsonts of - -r^ as measurements of (T,^ as opposed to measurements of C'f or J . 

da 1C 'L 1 

The HRR Field 

W<r nf^w give the HRR field in terms of (AT ) . Whereas similar relations 

have been written in terms of C* for ateady-state creep (36j, the relations in 

terms of (AT,) will be valid for nonsteady creep as well as steady-state creep. 

1 c 

The HRR field as given in (37] but modified for creep by replacing and u^ by 
and u^ respectively, As: 


^Note that these equations arc derived on the assumption that ^ij=0 (l.e., creep 
steady-state). Therefore, in order to have a well defined creep constitutive 
law we must have and n finite. 
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(T 


IJ 






zL 

n+l 


[5^j(0).8^^(e)i 


(2.49«) 


-n 

*1J " (2.49b) 

1 

u. - yK u,(e) (2,49c) 

1 £ 1 

where d (9) has been normalized to have a maximum value of unity and and 
eq ^ a 

are amplitude factors which are related by 

\ - (K^)" (2.49d) 

It can therefore be seen that the asymptotic crack-tip fields are entirely 
determined when (or K^) is known or specified. Combining (2.49a,b,c) with 
the first equality of (2.14), using (2.49d) and rearranging, gives: 


^0 — \ J \ yl* / 


1 

n+l 


(2.50) 


where Is analogous to 1 defined by Cq. (24) of [6] except for the factor 
n/(n+l) multiplying the energy density term. To be explicit, 


I* - I + 




(d (e)]""^^ cosOde 
eq 


(2.51) 


It is therefore seen that knowing the value of (Tj^)^ Is equivalent to knowing 
and thus is sufficient for defining all aspects of the asymptotic crack<tip 
field during nonstcady creep as well as under steady-state creep conditions. 


:i 
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SECTION III 


DERIVATION OF FINITE ELEMENT EQUATIONS 


Before etatlng the variational principle and deriving the finite elraent 
equationa, it ia beneficial to illuatrate the conatitutive behavior to which 
the finite element model ia addreaaed. 

Elaatlc/Viscoplaatlc Conatitutive Relatione 

A rather general rate conatitutive law propoaed by Perzyna [38], can be 
written in incremental form aa: 


Ac.- 


iJ - -L 




(3.1) 


where y, E and v are the elaatlc ahear roodulua. Young 'a modulua and Polaaon'a 
ratio, respectively, ^"^ij " ^ deviatorlc atreaa and y m 

viscosity constant of the material. In writing (3.1) it is implied that 


Ae 


IJ 


e vp 

Ac At. 


IJ 


where elastic and viscoplastic strain increments, res- 
pectively. The yield function governs the magnitude of the incremental 

ij «L 

vlscoplastlc strains through the function <^(F)> where 

0 for F ^ 0 
♦(F) for F > 0 


«t(F)> - 


(3.2) 
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The relative magnitudea for the incremental vlacoplaatic atraln eomponenta 
are aeen to depend on the factor Thla impliea that when con- 

sidered ae a vector in nine-dlmenaional atreas space, is always directed along 
the normal to the aurface 

By choosing F = f - o [a - (3/2(. . ' t ' such that 3f/^T ■ (3/2) 

eq «q ij IJ 13 

(x^j/o^^), and choosing ^(F) “ (F)" we find that (3.1) becomes 


Ae 


^ “ At' + Ax 6 +• (3/2 )y(o 

2m Ij 3E ^kk ij 


(3.3) 


This represents the special case of creep behavior which Is considered exclusively 
2 

in this study. It should be understood, however, that the flnice element 
model which is described below Is applicable to the more general behavior rep- 
resented by (3.1). 

Derivation of Finite Element Equations 

The finite element model is derived from the principle of virtual work 


f x,,6e,,dV - / t 6u.dS - 0 (3.4) 

J\j ^ Js ^ ^ 

a 

Ih the present finite element analysis, we assume only infinitesimal defor- 
mations and strains; hence there is no need to differentiate between the deformed 
and undeformed configurations. 

In writing (3.4) It should be noted that x^^ are the stresses 

existing at time t + At (where t is the current time), t^ are the prescribed 

tractions on S at t + At, and fiu.[6e - (l/2)(6u ,+6u. .)] are arbitrary 
o 1 ij i.j j*l 

compatible virtual displacements. 

Following customary procedures we introduce the element displacement 
shape functions which relate element displacement u^ to element nodal 


"'It is common to choose f=P in which case we have what is called an associative 
law. 

2 

While (3.3) is known to deviate from real material behavior (especially for primary 
creep) it is a widely used constitutive law and therefore has been adopted in the 
present study. 
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f 


rv 


- {u} - [Nj{q); 6u^ - {6u) • [Nj{6q) (3.5) 

We again use the customary notation wherein strain (and stress) components 
are placed In one-dimensional arrays 

{«} - fB){q}; {6«) - fB](6q} (3,6) 

Substituting (3.5) and (3.6) into (3.4) and applying conventional procedures 
for assembling element matrices into global matrices we have 


loLCX 


- J {t}'‘'lN]ds]|| Jg{«q}) 5 {F)^{6Q} - 0 
e 

Since {6q} are arbitrary virtual nodal displacements, it follows that 


lelef Jv - J* {t}^[N]ds]j = {F>^ - (0)^ (3.7) 

°e 

We now express the stresses {t} at t + 4t in terms of the current 
stresses, and the Incremental stresses corresponding to the time in- 

crement 4t: 

{t} 5 - {t}j + (3.8) 


In (3.8) and in the following, the I and I-fl subscrlptn designate the 
incremental solution with which the quantity is associated. Application of 
the incremental elastic constitutive law results in 


^^^I+l " ^^^el^I+1 

- (T), 


(3.9) 


where the incremental viscoplastic strains and [E] is the matrix 


”9 

• «« 


of elastic constants. Substituting (3.9) into (3.7), taking tho trsnspoas, 
and placing the known terns on the right hand side wa have the final form of 
the finite element equations: 

+ ^®vp^I+l ” 

where 


[K] - gig j [Bl'*^[E](BjdV 

'^e 

<">W ■ eL I 

“o 

e 

ts f [B]’’[E]{A€ }_.,dv 

vp I+l ele Jy vp I+l 


(3.11) 

(3.12) 


(3.13) 


. 

The above volume integrals are evaluated in the current work by 2x2 Gauss 
quadrature. Tie array {T} is input directly in terms of node point forces. 
Solution Procedures 

It should be noted that [K] of (3.10) is just the elastic stiffness 
and therefore only needs to be formed and decomposed^ once. This results in 
significant savings in the number of computations per time step as compared to 
methods using s*'iffness matrices which must be reformed at each step (i.e., 
tangent stiffness methods). It should also be noted that the term 
is computed from incremental viscoplastic strains which are esti- 
mated using and the material constitutive law (3.1). Only for the 

^The equations (3.10) are solved in the current work by the decomposition [K]"> 
where [D] is a diagonal matrix (the only nonzero entries are those 
on its diagonal) and [L] is a lower triangular matrix (the only nonzero entries 
are those below its diagonal); see for example [39]. 
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■pedal alcuation when the etreaaes do not change with time will this eatinate 
be exact. Having obtained the increfflental nodal dlaplacenents 
solving (3.10), one can easily find the total incremental straina via 

the incremental analogue of (3.6). We now describe two procedurea for ob> 
taining 

The first and simpler method to obtain i* substitute the esti- 
mated used in solving for (3.9). If one does this, thm 

it happens that 




and therefore (3.10) becomes for the next step*^ 


(3.15) 


[KHAQ}j^ 2 ■ ^®vp^I+2 ■ ^^^I+l 

- ^^vp^I+2 


(3.16) 


This method was compared to the following method and was found to require smaller 
time steps to achieve similar results. 

Rather than using the estimated values of (Ac and (3.19) to 

determine the constitutive relation (2.1) is Integrated over the cur- 

rent time step at each Gaussian quadrature point with the condition that total 
strain (e) varies linearly with respect to time from to (The 

present study uses an Eulerlan scheme with each time step being divided into 
five subincrements. ) The result of this procedure is better adherence to the 
postulated constitutive law at the expense of Introducing a somewhat unequlll- 
brated stress state. The amount of disequilibrium depends on the accuracy of 
the original estimate for the incremental viscoplastic strains and thus on the 
time step size. 

At this point one has two alternatives. The first Is to use the 


This procedure results In the current model reducing to that of Zlenklewicz 
and Corroeau [40]. 






vlscoplastlc scrain incranenta obtained through the time Integration procedure 
ea on improved estimate and to re-aolve (3.10) for the currant time atep. Thie 
procedure would, after aeverel iteretiona, raault in a atreaa atate which is 
equilibrated to within some small user specified tolerance. With this type of 
procedure the time steps could be as large as those used with tangent stiffness 
methods. Further, it is reasonable to expect the solution to be at least as 
accurate as if a tangent stiffness method were used.^ 

The second alternative ia to go Inmadiataly to the next time step with 
the understanding that the term {R}^ in (3.10) results in the disequilibrium 
from the Ith step being corrected in the 1+1 atep. This feature ia the result 
of the virtual work statement (3.4) being written in terms of total stress and 
tractions rather than incremental quantities. Owing to this corrective nature 
and to the diminishing returns one obtains from additional iterations, the 
second alternati/c is used in the present study. 

Regulation of Time Steps 

The creep calculations use a variable time step size which is auto- 
matically regulated by the finite element program based on two criteria. The 
first criterion is the maximum percent difference between the Incremental 
equivalent estimated creep strain and the Incremental equivalent integrated 
creep strain for all the Gauss points in the mesh: 


=• Max 


'^^EST ~ ^^INT 


Ae 


INT 


(3.17) 


The second criterion Is the maximum ratio of Incremental equivalent Integrated 


This procedure could actually be more accurate if similar constitutive law 
integration procedures and equilibrium iterations are not performed with the 
tangent stiffness procedure. Also, it has been shown [41] that many element 
types become overly stiff whsr using the tangent stiffness method for modeling 
const Itut ive behavior approaching incompressibility. This problem is not en- 
countered with the current method. 
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creep etrain Co Che uquivalenc elastic strain: 



The user specified, maximum permissible valuea for and are 
respectively. The size of the next step is then obtained from 


(3.18) 
and Cj, 


Ac 


I+l 


■ At, 



(3.19) 


Note that the initial time increment cannot be determined from (3.19) and must 
be specified by the user so as to satisfy the two step size criteria. 

In the present study, the values of and are 0.2 and 1.0, res- 
pectively. With these values, it has been found that the initial time steps 
are controlled by while later time steps are controlled by C^. The values 
of and are strongly affected by the mesh refinement since a finer mesh 
results in Gauss points being closer to the crack-tip and therefore having larger 
stresses and strain rates. To determine the sensitivity of the solution to the 
selection of and C^, a compact specimen was analyzed with the above criteria 
and also with and being halved (i.e., ■ 0.1 anu C2 ■ 0.5). It was 

found that the load point displacement differed by less than 0.5X for all 
time and that the steady-state solutions were essentially identical. It there- 
fore appears these values of and are small enough to ensure Chat the 
solutions to be discussed do not depend on these step size criteria. 


33 


SECTION IV 


ELEMENTS FOR SINGULAR CRACK-TIP BEHAVIOR 


This section describes and compares several two-dimensional crack-tip 
singularity elements. Perhaps the primary motivation for Introducing singu- 
lar crack-tip elements Into finite element modelu Is the significant savings in 

computational expense. It Is shown In Section V, for example, that 57 element 

- 1/2 

model with elastic, r , singular elements results In a more accurate solution 

than a non-alngular 102 element model. The savings In CP time in this case is 

greater than 50%. Generally, one must consider that some additional effort 

Is required to develop and implement a special crack-tip element and that this 

tends to offset the savings in CP time. It has been discovered, however, that 

the very commonly used, elght-noded. Isoparametric element can be made to 
- 1/2 

produce an r strain singularity by merely shifting mld-slde node locations 
via the node definition Input data [42, A3]. Therefore, a very convenient 
means for modeling linear clastic crack-tip behavior exists. It has also been 
shown [44] that a 1/r type strain singularity can he obtained with this element 
type thus providing a suitable element for non-hardening plasticity problems. 

For more general singularity behavior, such as the r strain singularity 

associated with the HRR crack-tip field of power-law plasticity or creep* 
one must resort to sepcially formulated elements. 

In the following sections, we consider special elements for linear 
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elastic problem (r strain singularity) as well as special elements for 
problems Involving the HRR, ^ strain singularity. In discussing 

these elements, an attempt is made to point out their advantages and dis- 
advantages. 

Elements for Linear Elastic Materials 
Although many special elements have been used for linear elastic frac- 
ture analysis (see Atlurl [45] for review), we consider here only the elght- 
noded isoparametric element. There are two basic forms In which the elght- 
noded Isoparametric element can be used as a linear elastic crack-tip element. 

In the first form, the two midside nodes adjacent to the corner node located at 
the crnck-cip are shifted toward the comer node so that they occupy the quarter- 
point of their respective sides. This form Is Illustrated In Fig. 4.1 by the 
T>po A crack-tip mesh. The second form in which the element con be used Is 
ilxustrated iit Fig. 4.1 by the Type B crack-tip mesh. In this form, the elght- 
noded element is degenerated to a triangular element by defining two comer 
nodes and their midside node to be the same node which is located at the crack 
tip . Then the two midside nodes adjacent to the crack-tip corner node are 
shifted to their quarter-points. It is important that only one node be used 
at the crack-tip, as opposed to three superposed nodes, since the latter case 
has been slwwn [44] to result in the 1/r type strain singularity. 

Barsouji [41,44] notes that numerical experimentation shows the degen- 
erate triangular form yields more accurate results than the nondegenerate 

element. He goes on to recommend that the four-sided configuration be abandoned 

- 1/2 

based on the premise that the r singularity exists only along the edges 
of the clement and also that the strain energy for this element becomes unbounded 
if exact integration is used [44]. Ying [46], on the other hand, concludes 
that the r singularity docs exist within the four-sided element as well as 
along its edges and that the strain energy for the element Is bounded (and thus 

- 3S - 
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Fi. 3 . 4.1 Several crack-tip mesh configurations 
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thti elenont atlffntaa la wall dafinad). 

Apart from tha abova conaidaratlona, thara ara two aapacta of tha 
doganaratc Type B elamant which Inharantly make It prafarabla to tha nondaganarata 
alaaant. The flrat la that tha procaaa of collapaing ona aida of tha alamant 
to coinclda with tha crack-tip raaulta In tha alamant local eoordlnataa balng 
tranafornad Into a form of polar eoordlnataa. Slnca the alaiaant'a ahapa 
functlona ara dafinad In tha alamant local eoordlnataa It la to ba axpactad 
that angular blaa will be much laaa apparent for thla element type. The aecond 
feature of the degenerate triangular element which makea 1 preferable la that 
It la geometrically better aulted for creating crack-tip finite element 
reeahea of arbitrary refinement. Since the angular dependence of the near tip 
aolutlon la aignlf leant, thla flexibility for Increaalng the meah refinement 
in the angular direction la Important. Combining these two aspecta with the 
numerical evidence cited by Barsoum, It aeema the triangular, degenerate 
element la the better element for modeling the elastic crack-tip singularity. 

For this reason, all quarter-point element calculations In the present study 
use the degenerate triangular form. In particular, the mesh configuration 
Type B of Fig. 4.1 has been used exclusively. 

Elements Cor Materials with HRR Crack-Tip Fields 

In the previous section, we diaCMgsed crack-tip elements for linear 

elastic material jehavlor. it. has been seen that the standard elght-noded 

- 1/2 

isoparametric element can be made to have the r strain singularity and thus 
Is useful for analysis of cracks in linear clastic materials. It can be shown 
(see for example, Atluri (45]) that this is the only singularity which this elght- 
noded clement can exhibit. However, it can also be shown that higher order 
elements of the Isoparametric family can result in singularities of the type 
r"n/(n+l)^ provided n is an integer. If we denote the order of the isoparametric 
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Intcrpolac Lun as m, than •inguUrltlaa of tha typa can ba cbtainad for 

e an intoger auch that tj<» {AS]. In tama of tha HRR povar law axponant n, 
thla maana we can hava alngularltiaa of tha typa whara n in an Intagar 

auch that £m-l. 

Baaad on tha abova dlacuaalon, It can ba aaan that it ia poaaibla to 
employ iaoparaaatic alemanta aa HRR crack-tip alananta providad ona ia aatiafiad 
with Intagar valuaa of tha powar law axponant » n. By chooaing tha highaat valua 
of n which one ia intaraatad in nodaling, ona can than prograai tha n4>l ordar 
iaoparamatric alamant. Tha problaa with thia approach ia that valuaa of n 
for common materiala can ba aa high aa 20. Thla Impliea that ona would nead to 
program an iaoparamatric eiemant of ordar 21. Whila thla ia parhapa within 
reason, it will be shown that non iaoparame trie elements can ba derived which are 
more readily implemented. 

Two Crack-Tip Eleracnta from tha Literature 

In this discussion of apacial elementa we limit consideration to two- 
dimensional, triangular elements with straight aides. Tha elements are derived 
in terms of the triangular polar coordinates (p,o) illustrated in Fig. A. 2 
and which arc related to the global cartesian coordinates (x^) by 

Xf - xj^ + p{x^ - xj) + j(o + l)(x^ - x|) } (A.l) 

In (A.l), the superscripts denote the node number. The crack-tip is assumed 
to be located at node 1 (i.e., at p*0) . The geometric mapping of (A.l) is 
similar to circular polar coordinates in that the transformation cannot be 
inverted at p-0. 

We now consider several choices for the assumed displacement fields 
within the triangular region. The first choice is 

u^(p,u) ■ + P^(u^ - u^) + j(o + l)(u^ - u^) } (A. 2) 
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-noded arbitrary strain singularity element 


«t 


It can be seen that (A. 2) la almllar In form to (4.1) except that p la replaced 
by P^(0<X<1) 80 aa to Induce a aingularlty in the dlaplacement derivatlvea. The 
three-noded element reaultlng from (4.1) and (4.2) la of thts type propoaed by 
Tracey and Cook [47]. Inapection of (4.2) reveala that thla element permlta 
rigid body tranalatlon but doea not permit rigid body rotationa or conatant 
strain modes. While problems exist for which thla element provides reasonable 
results, the lack of rigid body and conatant strain modes make this an un- 
desirable element for general anlaysls. 

We now consider a straightforward procedure which allows an alternative 
to (4.2) to be written which (i) provides all the rigid body modes as well aa 
all the constant strain modes, (11) results In the desired r^ displacement be- 
havior and (111) results In compatible displacement fields with adjacent 
elements. First note that we can be assured our assumed displacement field 
contains all rigid body and constant strain modes provided It can accomodate 
the following general displacement field 

Ui - ai + bj^Xj^ + Cj^x^ (4.3a) 

“2’“2*Vl*V2 

where a^, b^ and c^ are constants. Clearly, a^ and a^ provide for rigid 
translation modes, while b^ and c^ provide for constant strains, and ^ 22 * 

respectively. The constants Cj^ and b 2 provide a rigid rotation if Cj^ ■ -b 2 i 0 
and a constant shear strain, Cj^ 2 » ^1 * ^*2 substitute (4.1) 

into (4.3) and regroup terms we have the result 

- a* b*P + c*pa (4.4a) 

^^2 " ^2 ^ ^ 2 '^ ^ c*pa (4.4b) 

Starting from (4.4), we can now proceed to add terms as desired with the 
only condition being that we maintain compatible displacements with neighboring 
singular and nonsingular elements 
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We now coneldor a seneral approach for aatabllahlng dlaplacenuint ahapa 
functlona for triangular crack-tip alamanta. Since both dlaplacamant con- 
ponenta will follow the aame form, we drop the aubacrlpt for almpllclty. We 
now write 

u(o,o) - |(1 + o)fj(p) + i(l - o)f 2 <p) (4.5) 


with 


f2(P) ■ *2 

EjCp) ■ ®3 + bjP + CjP^ 


(4.6a) 

(4.6b) 


Inspection of (4.5) ahowa that on the element side l-2(u(p,-l) * ^ 2 ^^^ 
on aide l-3,u(p,l) ■ also be aeen that u(l,a) la linear on 

side 2-3. Since and f^ each have three unknowns, It follows that element 
sides 1-2 and 1-3 must have three nodes. This means two new nodes must be 
created. Since the geometric properties of the element do not depend on the 
locations of these nodes their positions along the edges of the element are 
arbitrary. In the following, however, we choose to place these nodes at the 
midsides. These new nodes correspond to positions 4 and 5 In Fig. 4.2. De- 
noting the nodal displacements by u** , J-1,5, we now use the following con- 
ditions to determine the six unknowns in (4.6a) and (4.6b) 


u(0,o) - u^; u(l,-l) ■ u^; u(I,l) ■ u^ 
u(j,-l) - u^; u(—,l) ■ u^ 


(4.7) 


The result Is 


1 u 2 

^2 " “3 " “ ’ ^2 " “ 


1 ,31 

u - c^; « u - u - c^ 


(4.8) 


2u'' - 
2^-' - I 


S ■ 


2u^ - 
2^-' - I 


Substituting (4.8) and (4.6) Into (4.5) and defining the functions multiplying 
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as wc have 

^ 5 1 

u(p,o) ■ I N.u (4.9) 

1-1 ^ 

where 


N, - 1 - p 


^ ^ m 5 . 5,5 

- ^(p - p) ; 


n: 




. 5,5 

V2' 


m5 .5,5 
N5 - ^2^2 


with 


0-2 


1-A 


1 


- j(l - a) + o) 

^ *^2 " 

By examination, It can be seen that (4.9) has terms similar to those of (4.4) 
and thus can represent all the rigid body and constant strain modes. From the 
form of (4. 5,6, 8) It can be seen that the element must be compatible with 
neighboring elements. Therefore, we have an element whlhc satisfies all the 
requirements which we originally stipulated. If the a in Fig. 4.2 and In (4.9) 
Is replaced by 2a* - 1 we recover Che form of the equations suggested by Stem 


[48]. 

As a result of the appearance of in the shape functions for this 
element, Che Integration Involved In evaluating the stiffness matrix (particularly 
the Integration with respect to P) Is not suited to Gauss quadrature. In [48], 
Stern derives a special Integration rule, which when combined with standard 
Gauss quadrature for Integration with respect to a, results In exact stiffness 
Integrations. Unfortunately, the expression (28) in [48] which specifies 
the relative radial location of Che quadrature points is apparently in error. 

The corrected expression is as follows: 


I 

X-1 


(4.10) 


1 + X - n 

2X - X(1 + X)X2”^J 


where and X 2 denote the values of p at which the sampling points are located. 
The corresponding weights are then given by 


X 


1 


1 X-1 

2 ^2 



w 


2 


1 

2 



^ 2*^1 


X-1 



(4.11) 


(4.12) 


X 2X»1 

This quadrature rule Integrates terms of the type P , P and P exactly. 

Since the rule has four parameters (x^,X 2 »w^»W 2 ) and is only required to Integ- 
rate three types of terms, the locations of the quadrature points are not 
uniquely defined. Selection of X 2 according to the following clrterlon results 
in both the numerator and denominator of the bracketed term of (4.10) being 
positive and thus results in a valid quadrature rule. 



(4.13) 


Stem [48] presents a family of elements which are developed so as to 
be compatible with surrounding nonsingular polynomial based elements of arbit- 
rary order. While it should be possible to verify that each member of this 
family does indeed satisfy the requirements which were discussed in deriving 
the above element, this procedure provides little Insight to the method for 
deriving such elements. In the next section, we generalize the procedure used 
in arriving at (4.9) to derive an element which is compatible with quadratic 
elements (e.g., eight-noded isoparametric elements). 

Derivation of a New Crack-Tip Element 

In this section we generalize the procedure used in the previous section 
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to derive a crack-tip clement which (1) contains all rigid body and constant 
strain modes, (ii) results in compatible displacement fields with neighboring 
singular and nonsingular elements and (iii) result? in an arbitrary strain singu- 
larity at the crack-tip fo the type r^”^ (0<X<1). 

We start in the same manner as before by writing the displacement field 
in the form 

u(p,o) - (l-o^)fj^(p) +^(a-l)f 2 (p) + ^(o+DfjCp) (4.14) 

where 

f^(p) - a^ + b^p + Cj^p^ (4,15) 

Noting that the form of (4.15) requires three displacements being specified 
along each of three radial line segments (o« -1, o • 0, a • 1) we Introduce 
node points at locations 4 through 7 as illustrated in Fig. 4.2. At this 
point it Is seen that this procedure will result in an interior node. It 
will be shown later that this node can be eliminated in a nuniber of ways. 

Denoting the nodal displacements by u , j"l,7 we now use the following 
conditions to determine the nine unknowns in (4.14,15): 

u(0,o) ■ u^; u(l,-l) - u^; u(l,l) - u^ (4.16) 

u(|-,-l) - u^; u(-|-,l) ■ u^; u(l,0) ■ u^ 

u(j,0) - u^ 


The result is 


“l " “2 “ ^3 " u 

6 1 2 1 
b^ - U - u - Cj^; b^ - U - u ~ C 2 ; b^ - u" 


(4.17) 


1 

u - c. 


1,, 7 6 1, 1/, 421. 

Cl - j(2u - u - u ); c^ - j(2u - u - u ) 

I/O 5 3 1. 

c^ - -(2u - u - u ) 
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with fi - 2* - 1 

Substituting (4.17) and (4.15) Into (4.14) «nd defining the functions Biultlplylng 
1 

u as we have 

^ 7 1 

u(p.o) - I nJu^ (4.18) 

i-1 ^ 


where 


with 


N? - 


, 1/ ^ % 

1 - p - -(p _ p) 

W ”3 " ^3^2* ”4 


,7,7 „7 11 „7 


♦2*1 




■ 1 - o ; 


2/ X . 
p(p - P) 


^(0 - 1); (j)j - “<J(o 1) 


Inspection of the derived shape functions verifies that this element satisfies 
all the requirements which we stipulated at the beginning of the derivation. It 
was noted that this element has an interior node. While interior nodes are 
generally avoided so as to reduce the bandwidth of the equations to be solved, 
it seems that in the case of crack-tip elements the advantage of having additional 
degrees of freedom in the vicinity of the crack-tip more than compensates for 
the few additional equations which are involved. 

We now consider several alternatives for eliminating the Interior node 

of this element and note that one of these results in the corresponding element 

of Stern (48]. We start by substituting (4,15) into (4.14); using the condition 

that f (0) ■ f_(0) • f-(0) and regrouping terms we can write 
1 2 J 

u(p,o) - Uj^ + + Cj^P^l (con’t on next page) 


(4.19) 






+ + (Cj - C2>P^1 



It should be understood that any arbitrary condition relating to one, several 

or all of through u^ will suffice to remove the interior node (i.e., node 7) 
from the element relations. However, it seams more natural to aliminata the 
node by the removal of one of the terms of (4.19). Noting that the first 
two terms of the first square bracket and the first term of the second square 
bracket represent the rigid body modes and constant strain modes, we are left 
with four terms which can possibly be deleted so as to eliminate the interior 
node. 

If we choose to eliminate the term c^p^ by constraining c^ to be zero 
we see from (4.17) that this implies ■ y(u^ + u^) . This choice would some- 
what defeat the objective of having singular displacement derivatives- and thus 
is not advisable. Furthermore, it Is inconsistent to retain the higher order 
terms op^ and o^p^ while not retaining p^. Note that the term op^ cannot be 

used to eliminate since u^ does not appear in its coefficient (l.e., in 

2 2 X 

either c^ or c^). Therefore we are left with the terms o p and o p . Either 

of these terms can be chosen to eliminate u'^. Stern's element [48] corres- 

2 

ponds to the case in which the coefficient of a p is set identically to zero. 

Of the elements discussed above, only the seven-noded element has been 
implemented in the present study. In Section V, this element is used for 
elastic analysis as well as for creep analysis. The special quadrature rule 
proposed by Stem [48] and summarized in (4.10) through (4.13) has been used * 
exclusively in evaluating this element's stiffness. 
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SECTION V 


CREEP CRACK GROWTH COMPUTATIONS 


Description of Problems 

The creep crack growth analyses which will be presented In this chapter 
deal with three distinct problems. The following s'ictlons Introduce each 
problem by describing the physical aspects such as geometry, loading 
and material properties as well as by describing why the problem was selected 
and what is hoped to be gained by Its consideration. 

All calculations In this chapter assume Infinitesimal strains and small 
deformations. The crack propagation calculations use quarter-point crack-tip 
elements and a mesh shift Ing/remeshlng procedure. 

Problem I; Non-Steady Creep of a Compact Specimen 

The compact specimen geometry wa<' chosen for study because of Its wide- 
spread use In fracture experiments and because numerical solutions for this 
problcn* have appeared in the literature thus providing results with which to 
compare. The dimensions of the specimen as well as the material properties 
and applied loading (see Fig. 5.1) were chosen to coincide with tlusse used 
recently by Ehlers and Riedel [33]. The problem Is used for a mesh refinement 
sensitivity study and for exploring various aspects of the C* con- 

tour Integrals during both nonstcady and steady-state creep. 

Several finite element meshes have been uased in the analysis. All 
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Fi". 5.1 Compact specimen geometry, loading 
condition, and material properties 
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of these meshes employ two-dimensional, sight -noded, IsoparsMtric elements. 

The Integrations for these elements are accoBq>llshed with 2x2 Gauss quadrsturs 
and therefore only elements with straight sides are employed. As seen from 
the meshes In Fig. S, 2, the pin-load Ing hole is not modeled. In all models 
the horizontal placement of the point load corresponds with the load line of the 
ASTM standard geometry (x • 25.0 mm). The vertical position is y ■ 32.5mm. 

A sensitivity study showed that shifting the load to y > 40 an has virtually 
no effect on the pertinent aspects of the solution. 

Most of the meshes contain collapsed elght-noded isoparametric elements 

at the crack-tip as Illustrated In Fig. 6.2. In several calculations, the 

midside nodes of these crack-tip elements are shifted to their quarter-points 

- 1/2 

so as to produce an r strain singularity at the crack-tip. Also, several 
calculations are performed with a special conforming seven-noded, *rlangular 
element which imposes the HRR, r type strain singularity. Tablfa 5.1 Iden- 

tifioB the meshes for which calculations are made and also gives Che load point 
displacement and J for the elastic solutions. These J values are compared Co 
those based on the expression for given by Srawley (49]. 

Problem II; Constant Velocity Propagation in a Creeping Strip 

This problem is concerned with a finite height, infinitely wide strip, 
with a semi-infinite crack. Loading consists of uniformly applied displacement 
rates at the top and bottom edges. This problem has been chosen for two 
reasons. First, since the strip is infinitely wide and the boundary conditions 
do not change with time, the propagating crack-tip fields can be expected to 
reach a "convecting steady-state'* creep condition. Here we use the phrase 
"convecting steady-state" to mean that the field remains unchanged in time with 
respect to a coordinate system which is centered at and moving with the crack- 
tip. This terminology is used so as not to confuse this condition with the usual 
steady-state creep condition in which material stress rates are zero. 
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a . The 55 (square crnck tip elements) ond 57 (triangular 

crock tip elements) element meshes 



b. The 102 element mesh (331 nodes; 642 do.f.) 

c. The 300 element mesh (941 nodes; l840d.o.f. ) 



5 .2 Finite element meshes for the compact specimen 
(contour integral paths are indicated by dashed lines) 
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In thu convocclng staAdy-atAt* caao, atrcaa ratea for MtarlAl ps>lnta 

arc not zero. An a reault, C* cannot (In a atrlct acmaa) be a valid crack-tip 

paramator. Aa ahould be clear fro« earlier dlacuaalona, (T.) la a valid 

1 c 

parameter at all crack apaada alnce It doea not require material atraaa ratea 
to be zero. From a practical point of view, If the crack apeed la low enough, 
then one can expect to be a uaeful quantity. By varying the crack propagation 
apeeda for thla problem over the range of velocitlea obaerved experimentally 
(fur a given material), It ahould be poaalbla to determine If thla range haa 
any portlona In common with the range of velocitlea for which C* la a uaeful 
parameter. 

The second reason for choosing this problem la that C* can be evaluated 
analytically for the special case when the crack la stationary. This allows an 
Im ependent check on the finite element calculations and serves aa a reference 
for the analyses In which the crack la propagating. The analytical evaluation 
of C* parallels the evaluation of for a similar elastic strip problsm as 
discussed by Rice [50]. (See Appendix E) It should be noted that C* has been 

t 

shown to be related to the steady-state value of (T.) and therefore It la 

1 c 

possible to obtain stationary crack cose from C* and equation 

(B.l) of Appendix B. The direct evaluation of (T.) in terms of either Its In- 

1 c 

tegral representation or Its energy representation requires knowledge of the 
stresses In the region of the strip adjacent to the crack-tip and therefore Is 
not a trivial task. 

The material properties used In this problem arc representative of 304 
stainless steel at 630^C. These material properties and the finite element 
discretization are given in Fig. 5.3. The mesh for chrs problem may at first 
appear rather coarse; however, elastic and steady-state creep solutions obtained 
with this mesh are sufficiently accurate to Justify its use. The comparison 
of computed clastic values and steady-state C* values with their analytic 
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values Is given in Table 5.2. 

Problem III; Creep Crack Growth in Double-Edgo-Crack Spsciraens 

The purpose of considering this problem is to apply the model to e 
problem for which experimental data exists. While much experimental date hes 
been reported in the literature, most authors do not Include sufficient in- 
formation to allow a numerical simulation of their experiments. The current 
problem is based on the experiments of Koterezawa and Iwata [51]. The primary 
reasons for selecting this work for study are that crack length versus time 
t.istories were given and that the experiments were performed on 30A stainless 
stool for which high temperature elastic and creep properties wore already 
available. 

The goomotry of the experimental specimens is given in Fig. 5.4. The 

finite element mesh for the calculations is shown in Fig. 5.5 with contour 

integral paths being indicated by dashed lines. It can be seen that the mesh 

takes advantage of the two planes of symmetry for the specimen md does not 
o X 

model the 60 notch. The initial crack length indicated in Fig. S.S corresponds 
to the notch depth in the spccime.. All calculations for this specimen assume 
plane stress conditions and use the material properties given in Fig. 5.3. 
Elastic results for two crack lengths are compared in Table 5.3 with those 
(based on formulas for K^) from [52] and are seen to be in good agreement. 

Compact Specimen Analyses 

The following describes several calculations for a compact specimen 
during transition from an Initial elastic state to one of steady-state creep. 

The geometry, loading, material properties and other details were described In 
the first section of this chapter as Probelm I. We first consider results for 

s 

the 300 element mesh of Fig. 5.2 in terms of C? and then (T.) . Then we address 

^Modeling the notch would have required the mesh shifting subroutines to be 
generalized. 
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(dimensions in mm) 


.4 Double-edge-notch specimen geometry used by 
Koteraiawa and Iwata [51]. 
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Fig. 5.5 Finite element mesh for simulation of 
doublo-edge-notch experiments 


(contour integral paths are indicated by dashed lines) 
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Table 5.3 Computational Aspects of the Elastic and Non -Steady Creep 
Portion of the Double -Edge -Crack Calculations 
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the topics of mesh refinement and the use of special crack-tip singularity ele- 
ments. 

Calculation of (t. ) and C* for Monsteady Creep 
1 c i 

The path- independence of (C|)^ during nonsteady creep is Illustrated 
In Fig. 5.6 using results from the 300 element mesh. The C superscript 
designates the particular contour which Is used, with 4 being the 
nondlmensional distance from the crack-tip to the point where the contour 
crosses the crack plane. Therefore, ^ Is zero at the crack-tip and has a 
maximum value of unity when the contour Is at the boundary of the specimen. 

r 

Values of (C^) are plotted as a function of time for nine values of C ranging 
from 0.03 to 0.92. It Is seen that (C*)^ is largest for contours close to the 
crack-tip (small O and that as steady-state Is approached, the values from 
all contours converge to C*. The solution has essentially reached steady- 
state at 300 hours. After 300 hours, the values of (C|)^ for all nine contours 
are within 1.5 percent of their average value. This value of C^, as well as 
values from calculations with the other meshes, Is given In Table 5.1. 

Now we consider computed values of approximated by (AT^)^/At. 

The values of (ATj^)^ are obtained through the specialization of (2.10) to the 
case of infinitesimal strains, small deformations, symmetric mode I behavior 
and traction-free surfaces: 


(«i)^ 


f [n^AW - „j(,y 
23A ^ 


. Aty) 



(5.1) 




3t 




Ae^dV 
3x, ij 


Since (T,)- is the limit of (f,)^ as c goes to zero, (T,)*^ is plotted as a 
function of e for several times (see Fig. 5.7). In this figure, e is the 
nondlmensional size of and is measured in the same manner as the 
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t = l0.8hr 





nondlmena tonal aize of ^ 234 * open points are the values of as co»- 

puted by (5.1), for nine contours in the 300 element model. The value of the 
crack-tip parimetor (T^)^ is given by the intersection of each respective curve 
with the e - 0 axis. Due to the large gradient in (t,)^ for small e it is 
seen that the accuracy of any extrapolation based solely on these evaluations 

s 

of (i.e., open points) would be of questionable accuracy, except perhaps 

near steady-state conditions. The solid points at c ■ 0 in Fig. 5.7 have been 

obtained using (2.14). It is seen that these values of (T,) appear to be reasonable 

1 c 

extrapolations of the curves of (T^)^ (5.1) thus giving some degree of confidence in their 
accuracy. 

Path Independence of (f.) 

1 c 

Based on arguments put forth in earlier portions of this paper, the 

value of (Tj^)^ obtained through (2.14) sliould be independent of the path 

(i.e., r^^^) which is used in its computation. This path-independence is 

illustrated by Fig. 3.8a where Is plotted as a function of the nondimen- 

sional distance of from the crack-tip, C, for several times. Generally, 

the path- independence is seen to be quite good. The largest deviation from 

path- independence in this figure is for the intermediate time of 10.8 hours 

with the difference between the extreme contour values being less than three 

percent. To further emphasize this path-independence, (Tj^)^ is plotted as a 

function of time in Fig. 5.6. As a result of its path-independence, (T.) is 

1 c 

represented by a single curve. Interestingly, this curve is a straight line 
for times before approximately 10 hours. 

Riedel and Rice [36] have arrived at the following approximation for 
K^(which they call A(t) based on the assumed approximate path-independence of 
during the initial portion of nonsteady creep: 


K 

0 


Kj(l - v^/e" 
(n+Dylt 


1 

n-fl 


(5.2) 
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b. 102 element mesh with 
quorter-point singularity 
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c. 57 element mesh with 
quorter-point singularity 


F15. 5.3 Path- independence of (T. ) , for several times 

during non steady creep 
is the non-dimensional size of 


63 - 


:as s; .*« 


Comparing (5.2) with (2.50) ono concludoi Chat ahould behava Ilka I/c for 

clmea when (5.2) la valid. In a Log«log plot of veraua tine thia would 

roault In a atralght line with a alope of -1. The atralght line ahown In Fig. 

5.6 1j Inclined from the horizontal by 40** and thar<^fore haa a alopa of -0.84. 

The current work haa reaultad In aome evidence that la approximately path- 
independent during Initial nonateady creep but that lea value tenda to Increaae 
with time. Thle tendency for to Increaee with time could explain the rather 
significant departure of the current reaulta from the behavior of (5.2). 
Quarter-Point Singularity Element Calculation 

We next consider the reaulta of computations uaing 57 and 102 element 
meshes with quarter-point slngularlclea. The purpose of considering these loss 
refined meshes Is to determine if the expense and effort In using the 300 element 
model Is necessary for obtaining accurate results. Table 5.1 summarizes Che 
results of these meshes for the limiting cases of purely elastic behavior 
and steady-state creep behavior. For the elastic problem it is seen t'.\at the 
results from these meshes agree with the 300 element mesh results co within 
one percent. At steady-state the 102 element model still agrees with the 300 
element mesh (in terms of Cp to within one percent while the 57 element 
model now differs by approximately eight percent. 

The contours used for the 57 and 102 element mesh are indicated in Fig. 
5.2. The 57 element mesh has four contours while the 102 element mesh has 
eight. The path- independence of computed from (2.14), is Illustrated 

for those two meshes in Fig. 5.8b and 5.8c. It is seen that the degree of 
path- independence In both Is similar to Chat observed for the 300 element mesh. 
Since we have evidence that die 57 element mesh is less accurate than Che other 
meshes at stcady'-state, it appears that high quality of Che path-independence 
cannot generally be Interpreted as meaning the solution Is accurate. 
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To dotermino th« ad«qu«cy of the 37 and 102 aleaent aaahaa for the non* 
ataady creep problea we now coapara their hiatorlea with that obtained 

with the .100 element meah (aee Fig. 3.9) . The curve appearing in thla figure 
haa bean placed through computed point a from the 300 element meah. The raaulta 
of the 102 element meah agree almoat perfectly with thla curve for tinea between 
0.2 houra and 16 houra. Prior to thla period and after thla period the reaulta 
fall below the curve by aa much aa 20 percent. While little can be aald about 
the abaolute accuracy of the calculationa for early portiona of nonateady creep, 

a 

we know (baaed on Appendix B) that ahould agree numerically with C* at 

ateady-atate to within a few percent. Therefore It can be aald that the 
valuea of from the 102 element meah arc algnlf Icantly In error at ateady* 

atate. Recalling that thla model gave a ateady-atate value of C* which agree 
quite well with the 300 element meah reaulta (aee Table 5.1) It la perhapa aur- 
prlalng that auch a algnlf leant error In the ateady-atate value of can 

exist. To better underatand the reaulta of thla model, plotted aa a 

function of in Fig. 5.10. It la noted from thla figure that the valuea of 
(T^)^ baaed on (2.1A) (l.e., the solid points) appear to be reasonable extrapo* 
latlons for times when the results are In agreement with the 300 element mesh 
results. However, as steady-state is approached, it is seen that these solid 
points no longer appear reasonable. If, however, one extrapolates the values 
of to e " 0 for the bottom two curves of Fig. S.IO, it is found that 

these values of are in good agreement with the 300 element mesh results. 

In comparing the equations for evaluating CJ, and IC la 

seen that la ihe only one of the three which involves an integration 

over the crack-tip quarter-point elements. Based on this and the apparently 
good accuracy of C* and (f,)^ It is believed that the integration over 
these elements is the major cause of discrepancy in between the 102 

and 300 clement mesh calculations. 
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Fig. 5.9 Comparison of (T, )_ from calculations 

X c 

with throe finite element meshes 
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t»0.309 hr 


t*4.27 hr 




t«26.lhr 


t*303 hr 


Fig. 5. 1C (T, as a function of € for several tines 

during non steady creep 

(results from 102 element mesh with quarter-point singularity) 







The 57 element reeulte do not compere fevoral 1y with the curve of Fig. 

5.9 for any eignlflcent portion of the solution. For most times the values of 

(T.) fell below the curve with the percent difference ranging from 50 percent 
X o 

mt t • 0.02 hr. to 15 percent et steedy-state. Based on the discrepancy of C| 
indicated in Table 5.1 and in the generally bad comparison of In Fig. 

5.9, it appears that the 57 element mesh with quarter-point singularity is not 
sufficiently refined for accurate creep calculations. This conclusion is per- 
haps a bit unexpected considering the degree of accuracy which this mesh displayed 
for the elastic probelm (see Table 5.1). Tlie reason for this drastic change 
of accuracy in going from elastic to creep behavior may be that the crack-tip 
strain singularity (l.e. , is inappropriate for the type be- 

havior expected to exist during creep. This topic is addressed in the following. 
HRR Singularity Element Calculations 

Based on the above observations, several analyses have been made using the 
seven-noded variable singularity element described previously (Section IV). The 
elastic solutions obtained with this element agree very well with those using 
the quarter-point Isoparametric element as can be seen from the entries of 
In Table 5.1. Also Included in Table S.l are the quasl-steady-state values 
of C*. The introduction of the correct strain singularity for steady-state 
creep (r does not significantly affect the 102 Element mode's Cjl but 

does Imporve that of the 57 element model. 

The analyses which use the seven-node singular element have the same 
singularity for the elastic solution and the subsequent creep solutions. At- 
tempts at changing the singularity from the elastic r to the r *'/^*'**’^^ 
value between the elastic and first creep solution have created numerical dif- 
ficulties due to the disequilibrium introduced in the process. No attempt 
at a gradual transition has been made. 


- 68 


It w 


original - 
OF POOR QUAinV 

The results using the seven-noded singular element are shown in 

Fig. 5.11. The solid curve represents the results of the 300 element model. 

The evaluation of (t.) is according to (2.1A) with the numerical procedures 

1 c 

being Identical to those employed with the quarter-point elements except for 
the contribution of the singular elements to the volume Integral. For the 
quarter-point elements, the stresses are assumed to be distributed llneraly with 
respect to the local coordinates. The volume integral is then evaluated in 
terms of quantities at the 2x2 Gauss ponts. For the seven-nude element, linear 
interpolation is used and in addition, several calculations are done assuming 
radial dependence of the type a - a + b seen from Fig. 5.11 

that none of the calculations agree well with the 300 element results. 


Based on this set of calculations, the general disagreement in 
between the singular crack-tip element models and the 300 element model does not 
appear to be due to the strength of the singularity which is introduced at 
the crack-tip. The general accuracy of C* for all the solutions with either 
elastic or creep type strain singularities supports this view. Rather, it seems 
likely that the difficulty in computing the volume Integral of (2.14) stems from 
the crack-tip element fields not satisfying the condition 


Lt Y -lillLlJ}! &e: (e,0)dQ - 0 


(5.3) 


From the discussion of Appendix A it can be seen that if this condition is not 
satisfied while at the same time the fields have the correct asymptotic (lingu- 
lar) radial dependence, then the volume integral of (2.14) does not exist. 

It therefore appears that accurate evaluation of (Tj^)^ using (2.14) 
cannot be accomplished if one uses crack-tip singular elements which provide 
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, ) as a function of time for several calculations with 
X o 

the variable singularity seven-noded, crack- tip element 


for the satisfaction of condition (S.3). At this point, it appears that ona must 

use a rather refined non-singular mesh (such as the current 300 element mesh) 

or introduce special crack-tip elements which satisfy (5.3) in order to compute 

(t^)^ accurately. The next section illustrates that for many probelms of 

practical Importance, a more attractive alternative may exist. 

(T.)^ as a Crack-Tip Field Parameter 
X c 

The previous discussion has pointed out some computational difficulties 
involved with evaluating vab concluded that these difficulties 

arc associated with the contlrbutlon of the crack-tip singularity elements to the 
volume integral of (5.1). It has been seen that despite the elastic strain 
singularity introduced by the quarter-point element scheme, the 102 element 
mesh gives accurate values of C*. Assuming this reflects the general accuracy 
of this solution. It Is desirable to use this relatively inexpensive model as 
opposed to using a very refined non-singular mesh oi !:o Introducing a special 
crack-tip element which satisfies condition (5.3). 

The effect of deleting the crack-tip singular elements from the volume 
Integral of (2.14) is shown in Fig. 5.12. Deleting these elements means that 
we are in fact evaluating where is the volume encompassed by the 

crack-tip elements. We will denote this particular as It will 

be shown that depending on the relative size of the crack-tip elements and the 

. 5 

proximity of the solution to steady-state condition, is a good approxi- 

mation to (T,) . 

i c 

The solid curve in Fig. 5.12 represents the results of the 300 element 

• 6 

mesh. The dashed curves are (T,) in the case of the 57 and 102 element 

1 c 

meshes and is (f,)^ with e « 0.03 in the case of the 300 element mesh. The 
i c 

errek-tip clement sizes for the 57 and 102 element meshes are 10 mm and 2.5 mm 
(or 20 and five percent of the ligament size) , respectively. 






• 6 

The values of (T.) for the 102 element mesh coincide with the solid 
1 c 

cuve for times after about 30 hours. Therefore !• • valid, path Indepen- 

dent, crack-tip parameter for times after about 30 hours and for values of (f.) 

1 c 

beginning at approximately 1.6 of the steady-state value. Fig. 5.6 shows that 
C* Is still significantly path- dependent at 30 hours and thus Is not an acceptable 
crack-tip parameter until much later. 

• g 

The curve of (T ) , (e ■ 0.03), for the 300 element model ^wems to Indl- 
1 c 

• 5 

cate that the validity of can be expanded to earlier times by reducing 

the size of the quarter-point elements. For example, a 5 of three percent of the 

. 5 

ligament would apparently result In being valid as early as seven hours and 

• 5 

for values of (T.) as large as 4,3 the steady-state value of (T.) . The curve of 
AC 1 c 

• (5 

(T^) ^ for the 57 element mesh tends to approach the solid curve as steady-state Is 

approached but never acutally converges even at steady-state. This indicates that 

• ^ 

this mesh Is too coarse for be a useful parameter. 

Constant Velocity Propagation In a Strip 
We now present some calculations for the cracked strip problem previously 
referred to as Problem XI. The geometry, loading and material properties for 
this problem are summarized In Fig. 5.3. The purpose of this problem Is to 
determine how significantly the crack-tip fields are affected by crack propa- 
gation velocities commonly observed in experiments. If for realistic crack 
speeds, the crack-tip field is essentially the same as for a stationary crack, 
then C| is path independent and characterizes the crack-tip fields. In any 
case, is a valid parameter. 

As noted previously, the steady-state C| values for the infinite strip 
problem can be obtained analytically without much difficulty. (See Appendix E.) 
Therefore, the procedure for this set of calculations is to select three values 
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of C* which span the range of valuaa reported In the literature. The valuaa 
which have bean choaan are 0.05. 5.0 and 50 N/mm.hr. Having thane valuea, the 
correaponding remote ateady>atate la determined aa well aa the edge dla- 
placenent which will reault in the tame remote elaatic Theee diaplacementa 

are applied to the model elaatically at t - 0. The reaulting valuaa of are 
compared to the analytic valuaa in Table 5.2. Next, the ateady-etate edge die- 
placement rates are determined. Using the elastic aolutlon as an initial 
state, the displacement rate, 8, is applied until the model reaches steady- 
state. The computed steady-state values of C* are compared to their analytic 
values in Table 5.2, The next step is to determine an upper bound crack 
velocity for each of the chosen values of C*. The following formula is based 
on the experimental data reported in [23,24] and represents data from center- 
crack, double-edge-crack, single-edge-crack, compact, and round-bar specimen 
types. 


da 

dt 


a C* 


1.173 


(5.4) 


where 

-2 

1.68 . 10 (upper bound) 
a - _ 

3.36 . 10 (average 

Having reached steady-state, the crack is propagated at the upper bound velocity 
given by (5.4) until it is determined that a convectlng steady-state has been 
reached. 


As noted previously, these calculations use the quarter-point crack-tip 
element. The crack growth simulation is accomplished through a combination of 
mesh shifting and remeshing as described in Appendix D. The nominal size of 
the crack growth increments is 0.4 mm or two percent of the crack-tip element 
width. For the highest velocity case (C* ■ 50 N/mm.hr), this procedure re- 
sults in crack growth at approximately every fifth solution step. 
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RoauXf for a Plane Strain Strip 


The ratults of th« plana atraln atrip calculation with C| • SO N/an.hr 

and ^ 1.65 nm/hr arc given In Fig. 5.13. Tha valuaa of *nd C* ara 

given for tha portion of tha calculation prior to ataady-atata aa wall aa 

during tha crack propagation protlon. Tha band rapreaanta tha range of valuaa 

obtained from tha four contoura llluatrated In Fig. 5.3. Both (T.) and C* 

1 c 1 

converge to the 50 N/mm.hr value at ateady-atate. During tha crack propagation, 

a 

it la aeen that and C* do not depart algnlf Icantly from their ateady- 

atate value. This means that this combination of loading and crack speed 

results in the crack-tip fields being essentially at steady-state conditions. 

» • 

This In turn means that both (T.) (or (T, ) ) and CF are valid crack-tip 

1 c 1 C 1 

field parameters. 

A closer view of the crack propagation portion of these curves Is given 
in Fig. 5.14. The dashed curves bracetlng the Initial portion of the solid 
curves represent the degree of path-independence and continue to be representative 
of the path- Independence observed during the crack propagation steps. For both 
(T^)^ and C*, it Is seen that the strip has essentially returned to Its steady- 
state condition prior to each crack growth increment. It is thought that the large 
departure of compared to C*) is more representative of the nonsteadiness of 

the crack-tip field, since the validity of C* in general, and the numerical evalua- 
tion of (2.43) in particular, are based on the existence of steady-state conditions. 
The effect of remeshlng Is seen at approximately eight hours. The 
first two steps after the remeshing were found to result In rather erratic con- 
tour integral values and are not indicated In these figures. The equilibrium 
correction feature of the present model and the automatic time step regulation 
procedure both act to quickly restore equilibrium at the crack-tip. 

The propagation portion of the calculation with C* ■ 5 N/mm.hr and 
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and C* for non steady creep, steady-state creep, and crack 
propagation phases of a plane strain strip analysis 




time, hr 

Fig. 5.14 CJ for the crack propagation portion 

of a plane strain strip analysis 
(Cj^ = 50 N/rr.m • hr ; v = I.65 T..T/hr) 
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- 0.111 mra/hr Is given In Fig. 5.15. Here sgsin it is seen that both (t,)^ 
at X c 

and C* have converged to the analytical value of C* (to within two percent, which 

is also about the degree of path->dependence) . Comporlng these i^esults with 

those in Fig. 5.14 for the higher C* and crack speed it is seen that steady- 

state creep conditions were not reached until 11 hours as opposed to approximately 

two hours in the previous case. Also, the return to the steady-state value after 

mesh shifting takes more time (two hours compared to 0.25 hours). However, when 

compared to the time between crack growth steps (both use 0.4 mm) it is seen 

that the lower velocity case return to steady-state well before the next growth 

step occurs. This result indicates that lower load levels and crack speeds 

are inherently closer to steady-state conditions. While this behavior may 

seem intuitively correct, it should be kept in mind that these results depend 

on the empirical formula (5.4) which is only valid for 304 stainless steel. 

It remains to be seen if similar beh. vlor occurs in other materials. 

A calculation has also been done for the case of C*">0.05 N/mm.hr. As 

a result of the large number of solution steps between crack growth steps, when 

-4 

using the maximum velocity of 5 . 10 mm/hr, the calculations used a higher 

_3 

velocity (5 . 10 mm/hr) . Even at this unrealisltcally high velocity (for this 

level of loading), the behavior is more steady-state-llke than the case of C*"5.0N/mm. 
hr described above. 

Results for a Plane Stress Strip 

In both plane strain problems discussed above, the steady-state value of 
(Tj^)^ is equal to C* to within the accuracy of the calculations. This is con- 
sistent with the relationship and comparison of C* and steady-state given 

in Appendix B. According to the approximate numerical values of this -<%ppendlx, 

• 

there should not be as close agreement between C* and (T.) in the case of 

1 i c 

plane stress. The primary purpose of this analysis is to verify this predicted 
behavior. 
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time, hr 

Fi;^. 5.15 (T, and C* for thy crack propagation portion 

of a plane strain strip analysis 
(CJ = 5.0 N/mm 'hr 5 v = 0 .111 mrr/hr) 







For this plana atraaa analyaia, C| waa choaan to ba 50 N/nm . hr and 
Che crack waa again propagatad at 1.65 an/hr. The remote t » the ateady-aCata 

yy 

uiaplacemcnt rate, 6, and Che elaatic diapiaceiMCPC , 6, ara 171 MPa, 0.168 nua/hr 
and 0.114 mm, reapecClvely. 

The reaulta of Chia calculation r.re given in Fig. 5.16 end 5.17. It 

• 6 

ia aeen from theae flgurea that (T.) doea converge to a aomewhat higher 

X c 

value at ateady-atate then C|. The ateady-ntate ia aeen from Fig. 5.17 to be 
approximately 52 N/naa . hr which ie higher chan by four percent. While 
chia ia a aomewhat amaller difference than auggeated by Appendix B, the aign of 
the difference ia the aame. In light of the approximate integration uaed in 
obtaining Che numeric valuea in the appendix, thia dlacrepancy ia within reaaon. 
As expected, Che general behavior for plane atreas conditions la essentially the 
same as for the previous plane strain analyses. Therefore, previous observations 
concerning the steady-state nature of the crack-tip field .during crack propa- 
gation are unchanged by Che shift to plane stress conditions. 

Double-Edge-Crack Specimen Analysis 

The following describes several calculations and their results for the 

problem previously referred to as Problem III. The geometry and finite element 

mesh for the double-cdgc-crack specimen are given in Figs. 5.4 and 5.5, respec- 

o 

tlvcly. The material properties are those of 304 stainless steel at 650 C and 
are assumed to be the same as those used in the strip analyses. (See Fig. 

5.3} Calculations have been made for remote applied stresses of 157 and 176 MPa. 
The experimental crack growth histories for these two stress levels are repro- 
duced from [51] In Fig. 5.18. It Is seen from these curves that the first two- 
thirds of the specimen lives arc characterized by crack velocities of less than 
0.01 mm/hr compared to nearly 0.5 mm/hr as rupture Is approached. 

The primary puipose of the following calculations Is to verify the 
conclusions which were reached In the previously described strip calculations: 
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and C* for non steady creep, steady-state creep, and 
crack propagation phases of a plane stress strip analysis 








that is, that the crack-tip field* are eaaontielly creep steady-state fields 
even Cor the most rapid creep crack velocities. These calculations will be a 
valid chock because the input to the calculations is only the remote applied stress 
and the measured crack velocity history, and does not in any way depend on experi- 

e 

mental determination of or as did the strip calculations. In fact, 

Kotcrazawa and Iwata do not report such measurements in (31]. 
analysis of Initial. Low Velocity Crack Growth 

This section describes the simulation of the initial portion of the 
crack velocity histories given in Fig. 5.18. In all of these calculations, the 
entire load is applied elastically at t - 0 and held constant throughout the 
subsequent creep solution steps. The convergence of and CJ to their 

steady-state values is sho\m in Fig. 5.19, with the dashed lines in the C* 
plots denoting the degree of path-dependence. It is seen that steady-state 
conditions are reached between a half and one hour after the load is applied. 

(Table ‘>.3 summarizes the computational aspects of this portion of the calcu- 
lation.) Therefore, it is seen by refereing to Fig. 5.18 that crack growth does 
not begin in the two specimens until well after steady-state conditions are 
reached. Since the current calculations assume small displacements and infinitesi- 
mal strains, and since only the strain and displacement magnitudes depend on time 
once steady-state is reached, there is no reason to continue the numerical calcu- 
lations to the crack initiation times indicated by the experimental results. 
Therefore, the initial crack propagation is simulated at times after steady-state 
conditions are reached but much earlier indicated by the experiments. 

The crack growth simulation results are shown in Fig. 5.19. The crack 
increment size for this study was approximately 0.01 mm which is nominally 2.4 
percent of the crack-tip element size. It can be seen that only one mesh shift 
(i.e., crack growth step) was modeled. It is clear from this figure that the 
time it takes for the specimen to return to steady-state is significantly less 
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than the time to the next crack growth Increment (Indicated by dashed lines). 
T)\erefore, the initial portion of the crack growth histories of Fig. 5.18 are 
clearly occurring under essentially steady-state conditions and thus C* as 

s 

well as valid crack-tip parameters. Since an Increase in C* re- 

sults in a more rapid return to steady-state conditions, the above conclusion 
will remain valid for the initial constant velocity portions of the curves of 
Fig. 5.18. 

When crack growth occurs so slowly that the crack-tip is essentially 
at steady-state, the crack-tip field does not depend on the history of the 
spcclmer. Therefore, assuming steady-state conditions continue to exist, it is 
possible to skip to the final stages of crack growth without modeling the 
Intermediate crack growth. If it Is found that crack growth is still slow 
enough for steady-state conditions to exist, then It seems reasonable to expect 
that the bevalvor at Intermelate crack lengths is also of a steady-state type. 

The following describes the results of this procedure when applied to the 
two double-edge-crack specimens. 

Analysis of Final Stage of Crack Growth 

To analyze the final stage of crack growth, the crack length is in- 
creased to 2.75 mm and the process of applying the load elastically and creeping 
to stcudy-statc is repeated. Table 5.3 summarizes the computational aspects of 
this process. The convergence of and C* to their steady-state values is 

shown in Fig. 5.20. Having reached steady-state, the cracks are gto*yn at the 
rate suggested by the. last portion of the crack histories (Fig. 5.18) as shown 

in Fig. 5.20, The significant increase in the frqucncy of mesh shifting (compared 

• 

to tViat in Fig. 5.19) due to the velocity increase makes the details of the 
curve difficult to distinguish in this figure. However, the time step size is 

such that six or seven steps occur between each crack growth increment. Unlike 

. 5 

the strip problejn, the values of (Tj^)^ and C* are clearly increasing during 
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. 5.20 ^1 analyses on final stage of 

crack propagation in two 
double -edge -crack specimens 


this crack propagation process. 

It Is necessary to determine whether this Increase In the crack-tip 
parameters Is due to the crack-tip no longer being at steady-state conditions 
or whether It Is due to the Increase In crack length. This Is accomplished 
by continuing the calculation without further crack extension. If the value 
of the parameters do not change significantly with time, this means the Increase 
was largely due to the crack length Increase and that crack growth Is still 
occurring under essentially steady-state conditions. Examination of the final 
portions of the curves of Fig. 5.20 shows that this Is the case. 

Based on this analysis; It appears that the conclusions reached as a 
result of the strip calculations are still valid. Since, (1) the strip analy- 
ses are much less expensive than this analysis of the double-edge-crack geometry, 
(ii) the steady-state C* for the strip Is easily obtained analytically and (ill) 
the crack-tip parameters do not depend on crack length for the strip geometry, 
it seems that similar studies for other materials and/or other temperatures could 
most effectively be accomplished through the use of the strip geometry. The 
need for such studies follows from the vast simplification of fracture analysis 
and prediction which results if crack growth ^ curs under steady-state conditions. 
More will be said about this point in the conclusions. 



SECTION VX 


CONCLUSIONS 


Summary of Result a 

A finite element model has been derived which Is generally applicable 
to viscoplastic material models. This model uses an Initial strain approach 
which reduces computation time spent on forming and decomposing stiffness 
matrices and also circumvents the problem of element Incompress lb llty constraints. 
Through special i'eatures, including a correction term in the finite element 
equation, this model provides for Improved adherence to the postulated constitu- 
tive behavior (as compared to the standard Initial strain approach) and 
allows time steps which approach in size those used in tangent stiffness 
methods. The accurocy and efficiency of this model with eight-node isopara- 
metric elements and the quarter-point crack-tip element approach have been 
verified through several calcu'’ations for a compact specimen geometry and a 
strip geometry. Also, a method of simulating crack growth through shifting 
of the quarter-point singularity elements and periodic remeshing has been 
described and demonstrated. 

It has been shown that despite the fact that C* characterizes the 
crack-tip fields under steady-state creep conditions, it does not have an 
energy or energy rate interpretation. A related path- independent Integral 
parameter however, does have the energy rate interpretation commonly 
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attributed to C^. Since experlmentalisti uia this energy Interpretation to 

e 

correlate creep crack growth rates. It appeara that (T.) (aa opposed to C*) 

X c X 

is gaining acceptance as a useful creep crack growth rate criterion . Further- 
more, does not rely on the existence of steady-state creep conditions 

and thus might be expected to be a valid criterion even if creep crack growth 
should occur at rates which preclude the existence of steady-state creep 
conditions at the crack-tip. 

A creep crack growth simulation for 304 stainless steel has shewn that 
for realistic load levels and corresponding crack speeds the crack-tip field 
Is essentially at a steady-state creep condition. This means that for this 
material, the propagating crack-tip field is largely unaffected by the 
history of crack growth or the history of loading. This feature can greatly 
reduce the analysis required for predicitng creep crack growth behavior 
in a component as can be seen from the following suggested methodology. 

clfii * 

We assume that t..c crack propagation speed — is related to (T.) 

• dt 1 css 

d £) 

(i.e., - -^) through the powe. law suggested by experimental data (23, 

24]. 


dt 1 css 


( 6 . 1 ) 


Next we determine (e.g., by steady-state creep finite element analysit^ 

« 

^^l^css * function of crack length. Because of the assumed steady-state 
crack-tip behavior, this can be accomplished by considering several discrete 
crack lengtlis and then fitting a curve. No crack growth simulation procedures 
are necessary. Combining (6.1) with this result provides the following 
relationship between time and crack length; 


■/. 


e(t) ((T.) ) 

1 css 


-0 


da + t 


( 6 . 2 ) 


where is the initial crack length and t^ is tlie time when crack growth initiates. 
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The only unknown quantity in (6.2) ia tha initiation tina t^. If tha initiation 
tlma for craap crack growth is aasumad to ba nagligibla (as might ba suggastad 
from the results of [9,10]) than (6.2) imnadiately providaa tha pradlctad crack 
growth history. 

Vitek 11] does not consider t^ to ba negligible based on several 
experiments (compact and doubla-edga-crack specimens) on two CrMoV steels. 

Using a dislocation model he further concludes that a measure of crack opening 
displacement (COD) correlates well with the initiation of crack growth in these 
experiments. If the same conclusion Is valid for 30''» stainless steel, then 
one can presumably predict t^ based on a transient finite element analysis of 
the initial flawed configuration and a critical value of COD. If Initiation 
occurs long after steady-state conditions are reached, it is then reasonable 
to estimate t^ using the rate of COD obtained from a steady-state finite 
clement solution. The use of (6.2) and of L.ie critical COD concept has not 
been investigated in this study. 

All of the creep calculations have used the constitutive law which is 
obtained by generalizing the Norton constitutive law to three dimensions. 

Whereas this law is a good representation of steady-state creep behavior, it 
does not, in general, represent the primary stage of creep. Future work should 
include a study of other creep constitutive laws (such as that of Bodner and 
Partom [32]). Also, the present model la derived on the assumption that dis- 
placements are small and strains Infinitesimal. The strains in the vicinity of 
the crack-tip for the present calculations with 304 stainless steel material 
properties arc on the order of 5-10!^ and therefore suggest that a finite 
strain formulation may be more appropriate. A study should be undertaken to 
examine this aspect of the model. 

As noted previously, the creep crack growth prediction methodology 
expressed in (6.1) and (6.2) has not been tested. A study to assess the 
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utility of this msthodology should thsrsfors bs undsrtsksn. This study 
should consider crack growth initiation as wall as crack propagation and 
should include a ranga of load lavals and savaral spaciman gaoMtrias. If 
the methodology is found to ba successful for constant applied loads, than 
tlie study should ba extended to consider more general load histories. 



APPENDIX A 


ExUttnca of Llaiti for Contour Int«ar>l Dtflnitlw 
This sppsndix dlacussss chs sxlstsnca of ths various Halts which hsvs 

s 

boon taksn In defining (dT) » (T) snd C*. In consldsrlni those Halts, we 

V V 

make use of the generally accepted result (see '6] for exaaple) that the 
strain energy density quantltlss W and W as well as the quantity W* behave 
as l/r In the vicinity of the crack-tip. This Is assumed to be valid for 
nonsteady as well as steady-state creep and also for the elastic state 
existing at t ■ 0. 

Baaed on the known asymptotic behavior at the eraek-tlp (l.e., the 
HRR fields) the limits of contour Integrals for equations (2.11,13,14,25,30, 
40,44) can be written In the following form provided one takes as being a 
circular contour centered at the crack-tip. 

J - f (A.L) 

J ^ 

A 

The nonsingular function f(c,9) becomes equal to f(0) when the limit Is taken 
and reflects the asymptotic nature of the HRR fields. 

In the following we limit the discussion to symmetric problems Involving 
only mode 1 crack-tip deformation. Further, we assume chat crack surface trac- 
tions and body forces ere Identically zero. With these conditions, we need only 

consider (AT,) and we can therefore rewrite (2.11) 

1 c 
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Of fOOR QUALITY 


^ f ln,4W - 
*'‘234 

. « r pi , 

«/« _M ®yi 

Now consider the limit of ths V^-V^ Intsgrsl of (A. 2). Inspection of this 
integral shows that It csn bs put into ths form: 


/v'-v ^ 


(A.?) 


t € 

R IT 


- C 




e)dedr 


where V* is s small circular region centered at the crack-tip and C is the 

la 

integral over the region - V^. The function ft(r»0) is a nonsingular 
function which becomes g(6) in the limit as r goeu to zero, where g(6) is 
known in terms of the HRR fields. Upon a first inspection of (A. 3) one is 
tempted to conclude that the limit does not exist since the Integrand has 
a non-integrable singularity Mt r - 0. If, however, we look at the right 
equality of (A. 2) it is seen that this conclusion results in a contradiction. 
Since we have shown that the limit of the integral on does exist (and 
therefore (A. 2) requires that the limit of the integral over must 

exist). A re-inspection of (A. 3) shows that the only way for this apparent 
contraction to be resolved is if the g(r,0) of (A. 3) has the following 

a 

property: 


n tt 

g(r,0)d0 g(0)d0 ■ 0 


(A.4). 


If function g(D) is luTown explicitly for the linear clastic case and there- 
fore (A.4) can be directly verified. For the URR field g(0) is not known 
explicitly and therefore (A.4) can only be verified numerically. 
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OF POOR QUALITY 


* 


For Inf InltoBlmal Mtrain, nonlln««r •Uitlclty, th« following rolation 
provides an oltomstlvo to vorlfylng (A. 4) dlroclty 

3u. 


V-V ^ 

€ 


L 


Jr 


n Lo — dS 
r J U 
‘234 




Ths rslstlon (A. 5) (which ssouaos csro crack aurfaco tractions and no body 
forest) illuatratas that this voluaa integral of type (A. 3) can be expraaaed 
in tarns of the contour integral of type (A.l). The relation (A.5) can be 
verified through the divergence theorem, the linear momentum balance condition 
and the following identities: 
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APPENDIX B 


Numerical Dlffaranca Betwaen (t, ) aad C* . 

A caa 1 

The purpoae of thia appendix ia to give aome examplea to illuatrate 
th«i numerical difference between given by (2.48). Uaing 

(2.50), (2.49a) and (2.51), we have 

- T - ‘ jr^5,^<6)=«ede (B.i) 

The valuea tabulated in Table B.I were computed approximately from valuee of 

I and plots of ( 8 ) given in [6] and should be viewed accordingly, 
eq 


Table B.I Comparison of (f.) and Cf 

X css X 

Plane Strain Plane Stress 

n"3n"13 n"3n«13 




0.98 1.00 


1.11 l.W 


It is seen that for the range of n commonly encountered, (t.) and Ct are 

1 css 1 

numerically very similar for plane strain but differ significantly for plane 
stress. 


C 
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APPENDIX C 


Numerical Mathoda for Evaluation of Cont>iur Infgralt 

The numerical procadurea for evaluating J. aa defined by (2.27), (dT.) 

X X c 

aa defined by (2.14) and Cj| aa defined by (2.31) are deacribad In the following. 
General Proccdurea 

In atudylng the contour Integral patha Indicated in the finite element 
meahea of Figs. 5. 2, 3, 5 (dashed lines) It la seen that the paths always pass 
through the centers of elements as opposed to along their edges. This procedure 
has been adopted so as to benefit from the presumably more accurate solution 
within the elements. Each element contour Is divided into two aegments with 
the integration being accomplished by two point Gaussian quadrature. All of 
the integrations arc performed In the element local coordinates. 

The J^-lntegral for Linear Elastic Analyses 

The contour integral portion of (2.27) Involves the stresses, , and 
the displacement derivatives, Both of these quantities can be evalu- 

ated at the required Gauss points through the element nodal displacements and 
simple manipulations with element matrices. In the current study, Is only 
considered as a parameter for linear clastic material behavior and therefore 


U - (l/2)cf^jC . 

The C*-Intcgral 


The C* Integral of (2.31) consists only of a contour integral. The W* 
I * 

of (2.31) Is evaluated using (2.43). The gradient rates are approximated by 

3Ui j DAu^ 

3x- ' At Dx. 

X 
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and therefore are average ratea for the increment as opposed to the rates at 
the end of the increment. The contour integration procedure for is as 

described above and uses two point quadrature for each element segment. Whereas 
stresses are easily computed at the required contour Gauss points in the elastic 
case, the stresses must be computed incrementally in creep analyses and there- 
fore stress information must be stored for each contour Integration point 
unless nonstandard element interpolations are used. In the present study, 
the St res . 08 at the contour Gau ^ points are Interpolated from the 2x2 element 
Gauss points through bilinear Lagranglan interpolation (in local coordinates), 
thus eliminating the need for additional storage. 

The ( Tj^) ^-Integral 

In the evaluation of (2.14) it is understood that are the stresses 
at the beginning of the time Increment being considered. The procedures for 
evaluating the contour integral portion of (2.14) are the same as used in evalu- 
ating C*. The incrmental stress-work density, AW, is computed from 


AW ■ ( T 


IJ 




The stress derivative appearing In the area integral of (2.14) is evaluated 
based on the 2x2 element Gauss point values and the assumption that the stresses 
are distributed blllnearly with respect to element local coordinates. Elements 
which are entirely within V^. are integrated with the usual 2x2 Gauss quadrature. 
Elements which are only partially within have each applicable quadrant in- 
tegrated by one point Gauss quadrature. 
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APPENDIX D 


Simulation of Crack Bxtenilon 

^toflollng 'ha propagation of a crack ualng the finite element method 
requires some special procedure for representing the creation of new crack 
surface. A common procedure Is to relax the nodal forces at the crack-tip 
node, thus In effect allowing the crack to extend to the next node along Its 
path of propagation. This relaxation process can be accomplished In one 
time step but usually la allowed to extend over several time steps due to the 
large change In nodal forces which Is inherent In the procesa. The major attrac- 
tion of this node-release procedure Is Its simplicity. There are two draw- 
backs of this procedure which resulted in an alternate procedure being adopted 
in this study. The first Is that the Increment in crack growth is directly 
determined by the nodal spacing In the mesh, therefore restricting the flex- 
ibility one has In selecting a time step size, the mesh size and /or the 
number of nodal force relaxation steps. The second and perhaps more important 
drawback Is that the method is not adaptable to models which use crack-tip 
singularity elements. 

A typical mesh in the vicinity of a crack-tip is shown In Fig. D.l. The 
region A represents the region being modeled by singular cisek-tip elements 
which in the present case remain centered on the crack-tip. The Type B 
elements are elght-noded isoparametric elciuents which distort so that the 
region A can remain centered on the crack-tip. The sequence of element configu- 
rations in Fig. D.l illustrates the shift ing/remcshlng procedure used In [53, 

5A] and adopted here. The region A is moved by shifting nodes without altering 
clement connectivity until the Type D element ahead of the region A becomes 
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Fig. D.l Examplo of mesh shif ting/remeshing procedure 
for simulation of crack growth 
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overly distoroil. At this point, the elemente in the vicinity arc redefined 
eo that further ahlftlng Is poaslble. It can be seen that this procedure allows 
the Increment In crack length to be arbitrarily small and does not Involve 
release of nodes In the same sense as for the previously described node- 
release procedure. 

The added flexibility «fforded by this shifting procedure does require 
some additional work. For example, In the creep crack growth application, nodal 
displacements and element Integration point stresses are Interpolated. The 
method of Interpolation which Is employed In this procedure Is discussed next. 

We consider that the solution at time t, has been obtained and we now 
must find the solution at time t^. During the Interval (b^,t 2 ) the crack has 
grown by an amount Aa. Since the crack growth simulation procedure requires 
that nodes be shifted, and since the solution at must be represented In terms 
of nodal and Gauss point quantities for the shifted mesh, it Is necessary to 
submit the affected nodes and Gauss point to an interpolation or fitting pro- 
cedure. 

The simplest Interpolation procedure for nodal displacements and the one 
used in [53,54] as well as for calculations in the present study Is one which 
directly uses the clement shape functions. In this method, the nodal positions 
for I e mesh at t^ are located in the mesh at t^^. Knowing which element 
of the mesh at t^^ encompasses this new node position allows the Immediate 
calculation of displacements by use of the element shape functions and the 
nodal quantities for the mesh at tj^. While this is a consistent procedure 
for transferring the solution at t^ to the mesh at t^, it should be understood 
that the transfer cannot be perfect. That this must be the case can be seen 
by considering that spatial derivatives of displacements, etc., are not continu- 
ous across element boundaries. Since the clement boundaries change position 
during the shifting process, points which had continuous derivatives at tj^ will 
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have dlicontlnous darlvatlvas In the meah at £2 vlco-veraa. 

When the meeh la ahiftad, tha alanent Gauaa polnta are alao ahlftad; thia 
means the Gauaa polnta rapraaant different aatarlal polnta before and aftar the 
ahlft. In order that the new Gauaa point atreaaae accurately rapraaant tha cur- 
rent stress state, it Is necassary to Interpolate stresses for tha new Gauaa 
point locations using the old Gauss point values and locations. The procedure 
for doing this la to assume the element stresses are distributed blllnearly with 
respect to the element local coordinates. Than It becomes possible to use 
bilinear Lagranglan Interpolation polynomials and the 2x2 element Gauss point 
stresses to Interpolate within each element. For all creep crack growth 
calculations in this study, the crick growth Increment sizes were chosen small 
enough that the new Gauss point stresses for each shifted element were always 
the result of interpolation within that same element. 


APPENDIX E 


Analytical Evluatlon of C* for th« Strip Probl— 


Thla app'^'^'llx briefly outllnaa the anlaytical evaluation of for two 
Infinite atrip probiema and then eumariaea the reealte In tabular fora. 

The flrat atrip problem la that which la llluatrated In Fig. 5.3. We ahall 
refer to thla probl«n aa Caa:!i A. The aecond problem, or Caae B. la alallar 
to Caae A in every reapect except the top and bottom edgea of the atrip are 
"clamped" rather than "on rollera". Theao boundary conditlona are aunoiarlaed 
aa follova: 

Caae A; u^(x,h) • -u^(x.-h) ■ t (E.l) 

Tjjy(x.h) - T^(x.-h) - 0 

Caae B: u^Cx.h) ■ -u^(x,-h) ■ S (E.2) 

u (x.h) - u (x.-h) - 0 

The crack surfacea are traction-free In both caaes. 

Wc can aelect a C*-lntegral contour which allows C* to be evaluated 
quite easily.^ Consider a contour of rectangular shape which coincides with the 
top and bottom edges of the strip, extends far enough ahead of the crack-tip 
so as to be in a steady-state stress field which is unaffected by the presence 
of the crack-tip, and extends far enough behind the crack-tip so aa to be In 
stress-free materia]. Wc now evaluate C*, as defined by (2.31), through the 
use of this contour. It can be sticn that for both Caie A and Case B, the 
horizontal portions of the contour at y ■ + h do not contribute to the 

^him procedure parallels that used by Rice [SOJ for the evaluation of J. in 
a similar clastic strip problem. 
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Integral, nor dooa the portion in the etreea-free naterlal. At the verltcel 
portion of the contour ahead of the crack-tip, the only non-sero term ia that 
involving W*. Therefore, it ia aeen that for both Caae A and Caae B we have 
CJ • 2W*h (E.3) 

where U* impllea W* exiating far ahead of tha crack-tip. Uaing the boundary 

ccnditlona (E.I) and (E.2) and the aaaumption of ateady-atate conditiona, it 

ia poaaible to evaluate the remote ateady-atate atreaaea, t^. Uaing (2.43) re- 

aulta in W* and* thua C*. The reaulta of thia exerciae are aunmarixed in Table 
• 1 

E.I. 

The correaponding linear elaatic atrip problem which ia obtained by re- 
placing the dlaplacement rate boundary conditiona by the correaponding dla- 
placement boundary onditiona hae been treated in a aimilar manner. Theae re- 
sulta are also given in Table E.I. 
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Table B .1 Analytical Solutlona for the 
Infinite Strip Problem 

Steady-State Creep Linear Elasticity 



Case A 


plane stress 11 1 

plane strain 
Case B 

plane stress 
plane strain 

This case does not have a steady-state solution since the 
boundary conditions require a volumetric strain rate. 
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